Ful lling Hard Constraints in Uncertain Linear ... - Semantic Scholar

Report 9 Downloads 21 Views
Ful lling Hard Constraints in Uncertain Linear Systems by Reference Managing Alberto Bemporad ? and Edoardo Mosca Dipartimento di Sistemi e Informatica, Universita di Firenze, Via di S. Marta 3, 50139 Firenze, Italy, tel.+39-55-4796258 fax +39-55-4796363,

fbemporad|[email protected]

Abstract A method based on conceptual tools of predictive control is described for tackling tracking problems of uncertain linear systems wherein pointwise-in-time input and/or state inequality constraints are present. The method consists of adding to a primal compensated system a nonlinear device called Predictive Reference Filter (PRF) which manipulates the desired reference in order to ful ll the prescribed constraints. Provided that an admissibility condition on the initial state is satis ed, the control scheme is proved to ful ll the constraints, as well as stability and set-point tracking requirements, for all systems whose impulse/step responses lie within given uncertainty ranges. Key words: Constraint satisfaction problems, Robustness, Nonlinear lters, Optimization problems, Reference input signals.

1 Introduction In applications the design of feedback controllers is often complicated by the presence of physical constraints: saturating actuators, temperatures and pres? Author to whom all correspondence should be addressed

sures within safety margins, working space limited by constructive restrictions, etc.. This issue has stimulated substantial theoretical advancements in the eld of feedback control of dynamic systems subject to input/state constraints (see Mayne and Polak (1993), and Sussmann et al. (1994), which also include relevant references, for an account of pertinent results). Most of this work has addressed the regulation problem in the presence of constraints |in particular input saturation| under the hypothesis that the plant model is exactly known. The main goal of the present paper is to address the constrained tracking problem for systems a ected by model uncertainties, by using tools from predictive control. Handling of hard constraints is in fact one of the potential bene ts of predictive control (Clarke, 1994; Keerthi and Gilbert, 1988; Mayne and Michalska, 1990; Mosca, 1995; Rawlings and Muske, 1993). Predictive control is based on the receding horizon control philosophy: a sequence a future control actions is chosen, by predicting the future evolution of the system, and applied to the plant until new measurements are available. Then, a new sequence is evaluated so as to replace the previous one. Each selected sequence is the result of an optimization procedure which takes into account two objectives: (i) maximize the tracking performance, and (ii) guarantee that the constraints are and will be ful lled|i.e., no \blind-alley" is entered. Recently, (Bemporad et al., 1997; Bemporad and Mosca, 1994) have applied receding horizon tools to the reference trajectory rather than to the control input, with a consequent substantial reduction of computational complexity. In fact, in contrast with other predictive control approaches (Kothare et al., 1996; Zheng and Morari, 1995), in (Bemporad et al., 1997; Bemporad and Mosca, 1994) the constraint ful llment problem is separated from stability, set-point tracking, and disturbance attenuation requirements, which in turn|in the absence of constraints|are supposed to be taken care of by a formerly designed compensator. By using the output maximal admissible sets theory (Gilbert et al., 1995; Gilbert and Tan, 1991) have developed a reference management technique for constraint ful lment in the ideal noiseless case. An extension to the case of input disturbances has appeared in (Gilbert and Kolmanovsky, 1995), where no model mismatch is considered. The aim of the present paper 1

is to lay down guidelines for synthesizing Predictive Reference Filters (PRF) for systems whose impulse and step response are uncertain in that they are only known to lie within given sets. The system is supposed to be a standard feedback loop, designed according to available robust control techniques so as to perform satisfactorily in the absence of constraints. Whenever necessary, the lter alters on line the input to the primal control system so as to avoid constraint violation and possibly maximize the tracking performance, according to a worst-case criterion. The paper is arranged as follows. Sect. 2 introduces the problem and presents the PRF algorithm for the class of systems under consideration. A description of the adopted model uncertainty is given in Sect. 3. Sect. 4 studies how to reduce the in nite number of constraints involved in the problem formulation into a nite number. Stability, tracking and other properties of the PRF are investigated in Sect. 5, while Sect. 6 is devoted to computational aspects. Finally, a simulative example of application of the PRF, which provides some design guidelines, is described in Sect. 7.

2 Predictive Reference Filter Design Consider a family S of linear asymptotically stable systems. Each member  of S has a state-space description of the form

:

8 > > > > x( > > > > > < > > > > > > > > > :

+ 1) = x( ) + Gg( )

y( ) = Hx( ) + Dg( )

(1)

c( ) = Hcx( ) + Dcg( );  2 S

where:  2 Z+ , f0; 1; : : :g; x( ) 2 R n is the state vector; g( ) 2 R p the command input, which in the absence of constraints would coincide with the desired output reference r( ); y( ) 2 R p the output which is required to track 2

r( ); and c( ) 2 R q the vector to be constrained within a given set C , which satis es the following

Assumption 1 C is a convex polyhedron with nonempty interior. Without loss of generality, we assume that C has the form

C = fc 2 R q : c  Bcg

(2)

In fact, a generic polyhedron described by inequalities of the form Acc  Bc can be rewritten in the form (2) by de ning a new vector c = Acc and, accordingly, new matrices Hc = AcHc, Dc = AcDc. Typically, (1) consists of an uncertain linear system under robustly stabilizing control. We take into account model uncertainties by assuming the true (unknown) plant is included in S , where S will be characterized in Sect. 3. Furthermore, inside the S , we choose a particular ^ called nominal system ^ :

8 > > > > ^( <x

^ ( ) + 1) = ^ x^( ) + Gg

> > > > :

c( ) = H^ cx^( ) + D^ cg( )

(3)

with x^ 2 R n . The aim of this paper is to design a Predictive Reference Filter (PRF), a device nalized to transform the desired reference r( ) to the command vector g( ) so as to possibly enforce the prescribed constraints c( ) 2 C ,  2 Z+, for all possible systems in S , and make the tracking error y( ) ? r( ) small. The ltering action operates in a predictive manner: at time  a virtual command sequence fg( ); g( + 1); : : :g is selected in such a way that, for all systems  2 S , the corresponding predicted c-evolution lies within C . Then, according to a receding horizon strategy, only the rst sample of the virtual sequence is applied at time  , a new virtual command sequence being recomputed at time  + 1. Several criteria (Bemporad et al., 1997; Bemporad and Mosca, 1994; Gilbert et al., 1995) can be used to select the class of virtual commands. 3

For reasons that will be clearer soon, we restrict our attention to the class of constant command sequences introduced in (Gilbert et al., 1995). This class is parameterized by the scalar , and each of its members de ned by

g(t +  j; ) , g( ? 1) + [r( ) ? g( ? 1)]; 8t 2 Z+:

(4)

At each time  , the free parameter is selected by the PRF via the optimization criterion

( ) ,

8 > > > > arg max > > > 2[0;1] 8 > > < > > > > < c(t +  ; ; x( ); ) > > > > subject to > > > > > > > > > : : t 0; 

j

2 C;

(5)

8  8 2S

where c(t +  j; ; x( ); ) denotes the predicted c-evolution at time t +  which results by applying the constant input g(t +  j; ) to  from state x( ) at time  onward. Then, according to the receding horizon strategy, at each time  the PRF selects

g( ) = g( j; ( )) Notice that requiring ( ) as close as possible to 1 corresponds to minimizing kg( ) ? r( )k, and consequently the norm of the tracking error ky( ) ? r( )k, depending on the tracking properties of (1). A scalar , or a constant command g 2 R p satisfying the constraints in (5) will be referred to as admissible at time .

Assumption 2 (Feasible Initial Condition) There exists a vector g(?1) 2 R p such that at time  = 0 the virtual command g (tj0; 0) = g (?1), 8t 2 Z+, is admissible.

For instance, Assumption 2 is satis ed for x(0) = (I ? )?1 Gg(?1), Hcx(0) + Dcg(?1) 2 C . Assumption 2, the particular structure of (4), and (5) ensure 4

that = 0 is admissible, and therefore the optimization problem (5) admits feasible solutions, at each time  2 Z+.

3 Model Uncertainty Description Uncertainty of dynamic systems models can be described in various ways. In the case at hand, frequency domain descriptions are not convenient because of the time-domain PRF design logic (4)-(5). Furthermore, if uncertainties involving state-space realizations are adopted, the e ect of matrix perturbations on the predicted evolutions become cumbersome to compute. Consider for instance a free response of the form (^ + ~ )tx(0): this gives rise to prediction perturbations which are nonlinear in the uncertain parameter ~ . On the contrary, uncertainties on the step-response or impulse-response samples provide a practical description in many applications, as they can be easily determined by experimental tests, and allow a reasonably simple way to compute predictions. Seemingly, step-response and impulse-response are equivalent, and one could be tempted to use either one or the other without distinction to describe model uncertainties. However, when used individually, both exhibit drawbacks. To show this, consider Fig. 1, which depicts perturbations expressed only in terms of the impulse response. The resulting step-response uncertainty turns out to be very large as t ! 1. However, this is not the case when each , for instance, contains an integrator in the feedback loop, which yields a unity DC-gain, and consequently vanishing step-response perturbations as t ! 1. Conversely, as depicted in Fig. 2, uncertainty expressed only in terms of the step response could lead to nonzero impulse-response samples at large values of t, for instance when the DC-gain from g to c is uncertain; therefore, any a priori information about asymptotic stability properties would be wasted. In order to minimize the conservatism of the approach (5), it is clear that the set  should be as small as possible. For this reason, in this paper we will jointly consider both step-response and impulse-response in order 5

1.5

Impulse response

1.5

1

1

0.5

0.5

0

0

0

10

20

0

Step response

10

20

Fig. 1. Step-response interval ranges (right) arising from an impulse-response description (left). 1.5

Impulse response

1.5

1

1

0.5

0.5

0

0

0

10

20

0

Step response

10

20

Fig. 2. Impulse-response interval ranges (left) arising from a step-response description (right).

to describe model uncertainty.

Assumption 3 Let  2 S and let H be the impulse response from g to c. Then, there exist a matrix M 2 R qp and a scalar , 0   < 1, such that, for all systems  2 S , jHtij j  M ij t; t 2 Z+; 8i = 1; : : : ; q; 8j = 1; : : : ; p

(6)

where Htij is the impulse response at time t from the j -th command input gj to the i-th constrained variable ci .

Notice that, although in (1) we are considering asymptotically stable systems, condition (6) characterizes only stability properties of the subspace which is observable from c. 6

The impulse response Ht can be expressed as the sum of a nominal impulse response 8 > > > ^ c^ t G^ > > > > :

if t > 0

D^ c if t = 0

and an additive perturbation H~ t. We describe the range intervals of H~ t as 8 > > > ~ tij > > > > :

2 [H ijt ; H ijt ] if t = 0; 1; : : : ; N ? 1

jH~ tij j  E ij t if t  N

(7)

where E 2 R qp , N is a xed integer, and i = 1; : : : ; q, j = 1; : : : ; p. In the same way, the step-response from g to c can be expressed as the sum of a nominal response

H is everywhere an underlined letter.

t?1

X W^ t , H^ c^ k G^ + D^ c

k=0

and an additive perturbation W~ t 8 > > > ~ tij > <W > > > > :

2 [W ijt ; W ijt ]

if t = 0; 1; : : : ; N ? 1

jW~ tij ? W~ tij?1 j  E ij t

if t  N:

(8)

4 Reduction to a Finite Number of Constraints Since the PRF operates over a semi-in nite prediction horizon, the optimization criterion in (5) involves an in nite number of constraints. In order to e ectively solve (5), we need to reduce this in nite number to a nite one. This will be achieved by borrowing techniques presented in (Gilbert and Tan, 1991), as follows. Under some assumptions on the desired reference r and the 7

W is everywhere an underlined letter.

past command inputs g(? ) applied before the PRF was switched on at time  = 0, Lemma 1 will show that the command sequence g( ) generated by the PRF is bounded. A new constraint on will be introduced, which ensures that in steady-state the predicted c-evolution stays away from the border of C . Then, Theorem 1 will prove the existence of a nite constraint horizon, the shortest prediction interval over which constraints must be checked in order to assert admissibility of a given virtual command sequence. Finally, an algorithm to nd such a constraint horizon will be provided.

Assumption 4 (Set-Point Conditioning) The reference signal r() satis es r( ) 2 R for all  2 Z+, where R is compact and convex. This amounts to assuming either that the class of references to be tracked is bounded, or that a clamping device is arti cially added to the PRF mechanism so as to satisfy Assumption 4. In practice, this is not a restriction since bounds on the reference are often dictated by the physical application.

Assumption 5 For all  > 0, g(? ) 2 R. Lemma 1 Provided that Assumptions 4 and 5 are satis ed, g( ) 2 R, 8 2

Z+.

r(0) g(-1)

g(0)

g(1) g(2)

r(1)

g(3)

r(2)

r(3)

Fig. 3. Reference set.

Proof. By (5), ( ) 2 [0; 1]. Then, as depicted in Fig 3, g( ) lies on the segment whose vertices are g( ? 1), r( ). By convexity of R, the result straightforwardly follows by induction. 2

8

In order to proceed further, we impose an additional constraint on the optimization (5). By (7), we can de ne W~ 1 , limt!1 W~ t, along with the relation

W~ 1 2 [W 1; W 1];

(9)

wherein W 1 , W N ?1 ? 1?N E , W 1 , W N ?1 + 1?N E , and, by Assumption 3, W^ 1 , limt!1 W^ t. Notice that (9) is not over-conservative in that there exist systems  2 S for which W~ 1 = W 1 or W 1. For an arbitrarily small  > 0, consider the following set n G = g 2 R : (W^ 1 + W~ 1)g  Bc ?  2 3

8W~ 1 2 [W 1; W 1]; 

67 6 7 6 7 6 7 = 666 ... 777 6 7 6 7 4 5



9 > > > > > > > > > =

2 R q >> > > > > > > > ;

(10)

For all constant command inputs g 2 G , the corresponding steady-state constrained vector cg , W1g is located in C at a distance from the border greater than or equal to a xed quantity, which depends on . The constraint g 2 G adds the following q additional constraints (W^ 1 + W~ 1)g  Bc ? ; 8W~ 1 2 [W 1; W 1]:

(11)

Hereafter, the constraints (11) will be added in (5) to determine ( ). In order to simplify the notation, in the next theorem we consider  = 0. Moreover, we replace x(0) with the sequence of past commands X? , fg(k)g1 k=?1, and de ne R , fX? : X?  Rg. Accordingly, c(t + 0j0; ; x(0); g) is denoted by c(t; ; X?; g), or, when the remaining arguments will be clear from the context, by c(t).

Theorem 1 Suppose that Assumptions 4 and 5 hold. Then there exists a nite 9

time t such that

c(t) 2 C ; 8t  t , c(t) 2 C ; 8t 2 Z+ 8 2 S ; 8X? 2 R; 8g 2 R;

(12)

2

Proof. See the Appendix.

Theorem 1 proves that in the optimization problem (5) it is sucient to take into account only the constraints up to time t . Since we are interested in the smallest time t , we de ne the following.

De nition 1 (Constraint Horizon) to = min ft : (12) holdsg t

(13)

Theorem 1 provides an existence result, as it proves that to in (13) is well de ned. However, in general, to is smaller than the quantity t estimated by Theorem 1 (see (32)-(33) for details). In order to nd out an algorithm which enables us to compute to , de ne Qt , f[; X?; g] 2 S  R  R : c(k; ; X?; g)  Bc, 8k = 0; : : : ; t, limj!1 c(j; ; X?; g)  Bc ? g. where, clearly, Qt+k  Qt .

Theorem 2 Suppose that Assumptions 4 and 5 hold. Then, Qt  Qt+1 ) Qt+k = Qt ; 8k  0

(14)

and

to = min ft : Qt = Qt+1 g t Proof. Let t such that (14) holds, and consider a generic [; X? ; g] 2 Qt+1 . By shifting X? = fg(?1); g(?2); : : :g in the new sequence X? , fg; g(?1); g(?2); : : :g, it follows that c(t + 1; ; X?; g) = c(t; ; X? ; g) Moreover, by Lemma 1,

10

X? 2 R . Then, [; X? ; g] 2 Qt , and being Qt  Qt+1 , [; X? ; g] 2 Qt+1 , which implies [; X? ; g] 2 Qt+2 . Hence, Qt+1  Qt+2 , or Qt+2 = Qt+1 = Qt . By induction, Qt+k = Qt , 8k  0. Let now tm , mintft : Qt = Qt+1 g. If [; X?; g] 2 Qto , then [; X?; g] 2 Qt , 8t  to, and hence to  tm . If, by contradiction, tm  to , then, by minimality of to, there exists [; X?; g] 2 Qtm 2 such that c(tm + 1; ; X?; g) 62 C , which implies Qtm 62 Qtm +1. Following an approach similar to the one used in (Bemporad et al., 1997) and (Gilbert and Tan, 1991), to can be determined by the following algorithm:

Algorithm 1 (Determination of to ) (1) t ?1 (2) Q?1 f[; X?; g] : c(?1; ; X?; g)  Bc, limj!1 c(j; ; X?; g)  Bc ? ,  2 S , X? 2 R , g 2 Rg (3) mit [;Xmax fci(t + 1; ; X?; g) ? Bci g, i = 1; : : : ; q ? ;g]2Qt (4) If mit  0, for all i = 1; : : : ; q, go to 7 (5) t t + 1 (6) Go to 3 (7) to t (8) Stop Observe that Algorithm 1 involves optimizations with respect to an in nite dimensional vector which contains X? and the impulse (or step) response coecients of system . However, by virtue of Assumption 3 (asymptotic stability), these can be approximated with arbitrary precision by nite dimensional optimizations. In fact, once a precision p has been xed, we can express the evolution of vector c as

c(t; ; X?; g) = Wtg +

M X k=t+1

Hk g(t ? k) +

with M such that 11

1 X k=M +1

Hk g(t ? k)

(15)

p 1 X X

j Hkij gi(t ? k)j  ip; k=M +1 j =1 8i = 1; : : : ; q; 8 2 S ; 8X? 2 R This allows one to implement Algorithm 1 by solving quadratically constrained quadratic programs (QCQP) with respect to [fg(?1); : : : ; g(?M )g; fH~ 0; : : : ; H~ M g; g]. Notice that M , and, consequently, the complexity of Algorithm 1, is related to the estimate  in (6).

5 Main Properties of PRF We investigate how the PRF a ects system stability when the reference to be tracked becomes constant. Since all systems  2 S are asymptotically stable, next Lemma 2 rst guarantees system stability by simply showing that g( ) converges to a vector g1. By using the viability result of Lemma 4, Lemma 5 will prove that, amongst the admissible command inputs, g1 is the vector which is closest to r. Lemma 6 will show that this limit is reached in a nite time. Finally, Theorem 3 will summarize the overall properties of the PRF.

Lemma 2 Suppose r( ) = r for all   0 . Then there exists g1 = lim !1 g ( ). Proof. If g(0 ) = r, it follows that ( ) = 1 is admissible for all   0 , and therefore g1 = r. Suppose g(0) 6= r. Since, by (4), g(0 + k) lies on a segment whose vertices are g(0) and r, by setting d( ) , kg( ) ? rk one has g( ) = r + kg(d0()?) rk [g(0) ? r] and 0  d( )  d( ? 1). Hence, since d( ) is monotonically non increasing and lowerbounded, there exists d1 , 2 lim !1 d( ), and, as a consequence, g1 = r + kg(d01)?rk [g(0) ? r]:

Lemma 3 Convexity of C implies convexity of the set G de ned in (10), 8 > 0. Proof. Consider  such that G 6= ; (G = ; is trivially convex), two set-points g0, g1 2 G , and g , g0 + [g1 ? g0], 0 < < 1. De ne W1 , W^ 1 + W~ 1.

12

Being C convex, the set C = fc 2 R q : c  Bc ? g is convex as well. Since W1g = W1g0 + [W1g1 ? W1g0], it follows that W1g 2 C , 8 2 S . 2 We show now that, given an admissible set-point g0 2 G and a new setpoint g1 2 G , there always exists a nite settling time after which a new admissible set point g is found by moving from g0 towards g1. In other words, as x( ) approaches the equilibrium state respective to g0, a new set-point in the direction of g1 becomes admissible.

Lemma 4 (Viability) Let g0, g1 2 G . At each time  there exist two positive reals and  such that the command g = g0 + (g1 ? g0) is admissible for all the past input sequences X? = fg( ? 1); g( ? 2); : : :g 2 R satisfying the condition

jgi( ? t) ? g0i j  ?t=2 ; 8t > 0; 8i = 1; : : : ; p

(16)

where g i denotes the i-th component of g.

2

Proof. See the Appendix.

Lemma 5 Suppose r( ) = r 2 R, 8  0. Then g1 = gr , where 8 > > > > < arg min gd gr = > d2R > > > : subject to gd

k ? r k2

, g(0 ? 1) + d[r ? g(0 ? 1)] 2 G :

(17)

Proof. Suppose by contradiction that g1 6= gr . We can apply Lemma 4 to the pair g1, gr . In fact, since R is bounded and g( ) 2 R, 8 2 Z, then there exist , 1 such that jgi(1 ? k) ? g1j  ?k=2, 8k  0, with  given by (36). By de ning as in (36), by Lemma 4 the command g = g1 + (gr ? g1) is admissible. Since d < d1 and d( ) is monotonically not increasing, d( )  d(1) = d < d1, 8  1, which contradicts d1 = lim !1 d( ). 2

Lemma 6 (Finite Stopping Time) If r( ) = r 2 R, 8  0, then there exists a nite stopping time s such that g( ) = gr , 8  s, with gr as in 13

(17). Proof. By Lemma 5, lim !1 g( ) = gr , and therefore Lemma 4 can be applied to the pair g(1), gr , with  given by (36) and 1 such that = 1 satis es (36).

2

Next Theorem 3 summarizes the properties of PRF.

Theorem 3 (PRF Properties) Suppose that Assumptions 1{5 hold, and that r( ) = r 2 R, 8  0 . Then, once the integer to is computed o -line via Algorithm 1, the optimization problem

8 > > > > arg max 2[0;1] > > > > 8 > > > > > > > > c(t +  ; ; x( ); ) > > < > > > > ( ) = > > < > > > subject to > t to ;  > > > > > > > > > > > > > > > > > > : : g ( 1) + [r ( ) g (

j

2 C;

8  8 2S ?

?

(18)

? 1)] 2 G

can be solved for all   0. By setting

g( ) = g( ? 1) + ( )[r( ) ? g( ? 1)]; 8  0; the constraints c( ) 2 C are ful lled for all   0, and for all  2 S . Moreover, after a nite stopping time s ,

g( ) = gr ; 8  s where gr is de ned in (17). In particular, if r 2 G , the PRF behaves as an all-pass lter for all   s.

Remark 1 Since after a nite time g( ) = gr , the asymptotical properties of the original system 1 remain unaltered, in particular x( ) ! (I ? )?1 Ggr as  ! 1. 14

6 Predictions and Computations To lighten the notation, assume  = 0, and consider the predicted evolution c(t) determined by a past command sequence X?  R and a future constant command g(t) = g, 8t  0, c(t) = P1 k=0 Hk g (t ? k) = c^(t) + c~(t), c^(t) , P1 ^ P1 ~ k=0 Hk g (t ? k), c~(t) , k=0 Hk g (t ? k). Equivalently, the nominal prediction c^(t) can be expressed in the computationally more preferable form,

c^(t) = H^ c^ t x^(0) + W^ t g;

(19)

for a consistent initial state x^(0). For the sake of simplicity, suppose that the quantity N in (7) is such that N > to , where to is the constraint horizon computed via Algorithm 1. Then, the prediction error c~(t) can be rewritten ~ k g(t ? k), 0  t  as c~(t) = W~ t g + PkN=?t1+1 H~ k g(t ? k) + c~p(t), c~p(t) , P1 k=N H to < N . We wish to obtain a recursive formula to determine the range of c~p(t) without requiring the storage of all past commands g(t ? k). Consider now a generic time  2 Z+ and de ne

cip(t +  j ) , =

max c~p(t +  j ) fH~ k 2 [?Ek ;Ek ]g1 k=N p X 1 X i

j =1 k=N

E ij k jgj (t ? k +  )j

At time  + 1, at the same prediction step t, one has

cip(t + ( + 1)j + 1) = cip(t +  j ) + p X E ij N gj (t j =1

j

? N +  + 1)j

(20)

Assuming that at time  = 0 system (1) is in an equilibrium condition, corresponding to g(?t)  g0 for all t > 0, (20) can be initialized with

cip(tj0) =

N

1?

15

p X ij g j E 0  j=1

j j

(21)

The constraint c(t+ j; ; x( ); g) 2 C can be nally expressed as the following constraint on g (W^ t + W~ t )g  Bc ? H^ c^ k x^( ) ? c~p(t +  j ) ? NX ?1 (W~ k ? W~ k?1)g(t ? k +  ) k=t+1

(22)

or equivalently, by (4), as the following constraint on (at + yt)  bt + xt ; t = 0; : : : ; to

(23)

where at , W^ t[r( ) ? g( ? 1)], bt , Bc ? H^ c^ t x^( ) ? W^ tg( ? 1), xt , c~p(t +  j )?PkN=?t1+2 (W~ k ?W~ k?1)g(t?k+ )?W~ t+1 g( ?1), and yt , W~ t [r( )?g( ?1)]. In the same manner, constraints (11) can be rearranged in the form (23). While vectors at and bt are known, vectors xt and yt are linear functions of the uncertainties W~ k . Observe that the particular form (4) has allowed to make xt , yt independent; in fact, yt is a function of W~ t, while xt depends on fW~ k gkN=?t1+1. Therefore, xit 2 [xit ; xit], yti 2 [yit ; yit], for each component i = 1; : : : ; q, and for all t = 0; : : : ; to, with 8
> > ~ kij > <W

> > > > ~ kij :W

2 [W ijk ; W ijk ] ? W~ kij?1

2 [H ijk ; H ijk ]

(24)

9 =

yit = min : W~ kij [rj (0) ? gj (t ? k)]; j =1

subject to W~ kij 2 [W ijk ; W ijk ]

(25)

and xit , yit de ned analogously. 16

x is everywhere an underlined letter. y is everywhere an underlined letter.

6.1 -Parameter Selection

The constraints involved in the optimization problem (18) can be rewritten as (ait + yti)  bit + xit ; 8xit 2 [xit; xit ]; 8yti 2 [yit ; yit ] (26)

t = 0; : : : ; to + 1; i = 1; : : : ; q where xit , xit, yit , yit are obtained by (24)-(25) for t  to , and the index t = to +1 has been de ned for the constraints deriving by (11). Consider the generic constraint (a + y)  b + x

(27)

By de ning

f (x; y) , ab ++ xy ; x 2 [x; x]; y 2 [y; y] n fag and recalling the result in Bazaraa and Shetty (1979), pp. 101, 473, f (x; y) assumes global minima and maxima on the extreme points of its domain. Consequently, it is easy to show that (23) has the following solution 8 > > > > > > > > > < b+x > a > > +y > > > > > > :

 minf ab++xy ; ab++xy g if y > ?a   ab++xy

(28)

if y < ?a < y

 minf ab++xy ; ab++xy g if y < ?a

Constraints (26) can be summarized as max

i = 1; : : : ; q

it  

min

i = 1; : : : ; q t = 0; : : : ; to + 1

t = 0; : : : ; to + 1

17

it

R V

θM βM

θL

T

JM ρ

JL βL

Fig. 4. Servomechanism model.

where, possibly, it = ?1, it = +1.

Remark 2. Recalling (24), the bounds xit , xit are evaluated as solutions of

linear programs. The computational burden can be hugely lightened if uncertainty ranges are given only on W~ t , t = 0; : : : ; N ? 1. In this case, (24) can be solved as trivially as (25) without the need of linear programs, at the cost of a more conservative description of system uncertainty, as pointed out in Sect. 3.

Remark 3. In the present formulation, the state x^( ) is calculated in an open-

loop manner by iterating (3). Apparently, this would lead to high sensitivity w.r.t. sensor noise, which in this paper is not taken into account. In fact, no feedback from  takes part in the selection of g( ). However, one should remind that system , in general, represents a feedback linear loop. For example, for a constant reference trajectory r( )  r 2 G , after a nite stopping time s system  receives g( )  r. Then, y( ) will track r with the error rejection properties deriving from the primal linear controller.

7 An Example The PRF is applied in connection with the position servomechanism schematically described in Fig. 4. This consists of a DC-motor, a gear-box, an elastic shaft and an uncertain load. Technical speci cations involve bounds on the 18

is everywhere an underlined letter

shaft torsional torque T as well as on the input voltage V . System parameters are reported in Table 1. Model uncertainties originate from the moment Table 1 Model parameters Symbol Value (MKS)

Meaning

LS

1.0

shaft length

dS

0.02

shaft diameter

JS

negligible

shaft inertia

JM

0.5

motor inertia

M

0.1

motor viscous friction coecient

R

20

resistance of armature

KT

10

motor constant



20

gear ratio

k J^L

1280.2

torsional rigidity

20JM

nominal load inertia

L

25

load viscous friction coecient

Ts

0.1

sampling time

of inertia JL of the load, which is determined by the speci c task. Denoting with M , L respectively the motor and the load angle, and by setting xp , [L _L M _M ]0 , the model can be described by the following state-space form 2

3

6 0 6 6 6 6 k 6 6 x_ p = 66 JL 6 6 0 6 6 6 4 k 

1

? ? JLL

"

JM

k JL

0

0

0 ? 2kJM ? #

L = 1 0 0 0 xp "

0

3

2

0

7 6 0 7 6 7 6 7 6 7 6 7 6 0 0 7 6 7 xp + 6 7 6 7 6 7 6 0 1 7 6 7 6 7 6 4 kT M +kT2 =R 5

JM

RJM

#

T = k 0 ? k 0 xp 19

7 7 7 7 7 7 7 7V 7 7 7 7 7 7 5

Since the steel shaft has nite shear strength, determined by a maximum admissible adm = 50N=mm2 , the torsional torque T must satisfy the constraint

jT j  78:5398 Nm

(29)

Moreover, the input DC voltage V has to be constrained within the range

jV j  220 V

(30)

The model is transformed in discrete time by sampling every Ts = 0:1s and using a zero-order holder on the input voltage. A robust digital controller is designed by pole-placement techniques, and has the following transfer function from e = (r ? L ) to V 9:7929z 3 ? 2:1860z 2 ? 7:2663z + 2:5556 ; Gc (z ) = 1000 4 (31) 10z ? 2:7282z 3 ? 3:5585z 2 ? 1:3029z ? 0:0853 The resulting closed-loop system exhibits a very fast response but inadmissible voltage inputs and torsional torques for the references of interest, as shown in Fig. 5 for a set-point r = 30 deg. The PRF is applied in order to ful ll (29), (30) for the uncertainty range which derives by an unknown load JL , 10JM  JL  30JM (the corresponding impulse/step response uncertainty limits were obtained by maximizing Htij ; Wtij with respect to JL). Fig. 6 shows the resulting uncertainty set for both the impulse and step responses. The constraint horizon is to = 15. A nominal load inertia J^L = 20JM is selected, along with N = 17, E = 1000[1 1 1 1]0 ,  = 0:8,  = 10?6. As a design rule of thumb, in order to have a description of the family of plant S as less conservative as possible, N should be approximately equal to the \length" of the impulse response in terms of time-steps: suciently large to describe accurately the range of variation of each sample when this is perceptibly nonzero, but also small enough to minimize computational complexity. Fig. 7 shows the resulting trajectories for a set-point r = 30 deg and a load JL = 25JM . This were obtained in 112s by using Matlab 4.2 on a 486 DX2/66 personal computer, with no particular care of 20

Load Position (deg)

60 40 20 0 0

0.5

1

1.5 Torsional Torque (Nm) Time (s) 200 1000 100

2

2.5

3

Input Voltage (V)

500

0 0

-100 -200 0

1 2 Time (s)

-500 0

3

1 2 Time (s)

3

Fig. 5. Unconstrained linear response. The shadowed area represents the admissible range.

code optimization. The standard Matlab LP.M routine was used to solve linear programs. Fig. 8 describes the e ect of the width of the uncertainty interval. The larger the uncertainty range, the more conservative the PRF action, and hence the slower the output response. In order to make comparisons, constraint ful lment is also achieved by linear control. The gain of controller (31) is reduced by a factor 16:9802 in order to have a maximum admissible set-point of 180 deg for the nominal plant J^L = 20JM (with such a gain the linear loop reaches the maximum admissible torque T during the transient). The resulting trajectories are depicted in Fig. 9 (thin line). In the same gure, the trajectories produced by applying the PRF together with the fast controller (31) are shown (thick lines). Fig. 9 shows also the resulting trajectories for r = 90 deg. While for the linear closed loop the rise-time is the same for both r = 180 deg and r = 90 deg (thin lines), this is no more true for responses obtained by applying the nonlinear PRF 21

Impulse Response (Torque) 400 200

Impulse Response (Voltage) 4000 2000

0 0

-200 -400 0 10 20 Step Response (Torque) 500

-2000 0 10 20 Step Response (Voltage) 3000 2000

0

1000 0

-500 0

10 Discrete time

-1000 0

20

10 Discrete time

20

Fig. 6. Uncertainty ranges for 10JM  JL  30JM (thick lines) and nominal J^L = 20JM response (thin lines).

(thick lines). In general, even if a linear controller gives the same performance of the PRF for the maximum admissible set-point, when the desired reference sequence is nonconstant a better tracking is provided by using a nonlinear reference lter.

Remark 5. As a general rule of thumb to design controllers which will be used

in connection with a PRF, in order to maximize the properties of tracking one should select a robust controller which provides fast closed-loop response for all the systems of the considered family. This usually corresponds to large violations of the constraints, which therefore can be enforced by inserting a PRF. On the other hand, this cannot improve poor tracking properties of the original system because, as observed in Remark 1, the PRF becomes an all-pass lter when the constraints are inactive. 22

Load Position and Generated Reference (deg) 40 30 20

θL(τ) g(τ)

10 0

0

0.5

1

1.5 Time (s)

Torsional Torque (Nm) 100

400

50

200

0

0

-50

-200

-100

0

1 2 Time (s)

-400

3

2

2.5

3

Input Voltage (V)

0

1 2 Time (s)

3

Fig. 7. Response for J^L = 20JM , 10JM  JL  30JM , and a real JL = 25JM .

8 Conclusions This paper has addressed the robust PRF problem, viz. the one of ltering the desired reference trajectory in such a way that an uncertain primal compensated control system can operate in a stable way with satisfactory tracking performance and no constraint violation in the face of plant impulse/step responses uncertainties. The computational burden turns out to be moderate because of the underlying simple constrained optimization problem.

Acknowledgments The authors gratefully acknowledge the constructive criticisms of the anonymous reviewers on the rst submitted version of this paper. 23

Load Positions (deg) 60 40 20 0 0

0.5

1

Torsional Torque (Nm)

1.5 Time (s)

100

400

50

200

0

0

-50

-200

-100 0

1 2 Time (s)

-400 0

3

2

2.5

3

Input Voltage (V)

1 2 Time (s)

3

Fig. 8. Response for J^L = JL = 20JM , and di erent uncertainty ranges: no uncertainty (thick solid line), [15JM ; 25JM ] (dashed line), and [2JM ; 40JM ] (thin solid line).

References Bazaraa, M.S. and C.M. Shetty (1979). Nonlinear Programming | Theory and Algorithms, John Wiley & Sons. Bemporad, A., Casavola, A. and E. Mosca (1997). Nonlinear control of constrained linear systems via predictive reference management. IEEE Trans. Automat. Control, AC-42, 340{349. Bemporad, A. and E. Mosca (1994). Constraint ful lment in feedback control via predictive reference management. Proc. 3rd IEEE Conf. on Control Applications, 1909{1914. 24

Load Position (deg)

200 150 100 50 0

0

100

2

4

6 Time (s)

Torsional Torque(Nm)

300

10

12

Input Voltage (V)

200

50

100

0

0

-50 -100 0

8

-100 2 4 Time (s)

-200 0

6

2

4 Time (s)

6

Fig. 9. Set-point r = 90; 180 deg, JL = 20JM , no uncertainty. Fast controller + PRF (thick lines) linear controller (thin lines).

Clarke, D. (1994). Advances in model-based predictive control, Advances in Model-Based Predictive Control, pp. 3{21, Oxford University Press Inc., New York. Gilbert, E. and I. Kolmanovsky (1995). Discrete-time reference governors for systems with state and control constraints and disturbance inputs. Proc. 34th IEEE Conf. on Decision and Control, 1189{1194. Gilbert, E., Kolmanovsky, I. and K. T. Tan (1995). Discrete-time reference governors and the nonlinear control of systems with state and control constraints. Int. Journal of Robust and Nonlinear Control, 487{504. Gilbert, E. and K. T. Tan (1991). Linear systems with state and control constraints: the theory and applications of maximal output admissible sets. IEEE Trans. Automat. Control, AC-36, 1008{1020. 25

Keerthi, S. and E. Gilbert (1988). Optimal in nite-horizon feedback control laws for a general class of constrained discrete-time systems: stability and moving-horizon approximations. J. Opt. Theory and Applications 57, 265{293. Kothare, M., Balakrishnan, V. and M. Morari (1996). Robust constrained model predictive control using linear matrix inequalities. Automatica, 32, 1361{1379. Mayne, D. and H. Michalska (1990). Receding horizon control of nonlinear systems. IEEE Trans. Automat. Control, AC-35, 814{824. Mayne, D. and E. Polak (1993). Optimization based design and control, Preprints 12th IFAC World Congress, 3, 129{138. Mosca, E. (1995). Optimal, Predictive, and Adaptive Control, Prentice Hall, Englewood Cli s, New York. Rawlings, J. and K. Muske (1993). The stability of constrained recedinghorizon control. IEEE Trans. Automat. Control, AC-38, 1512{1516. Sussmann, H., Sontag, E. and Y. Yang (1994). A general result on the stabilization of linear systems using bounded controls. IEEE Trans. Automat. Control, AC-39, 2411{2424. Zheng, A. and M. Morari (1995). Stability of model predictive control with mixed constraints. IEEE Trans. Automat. Control, AC-40, 1818{1823.

Appendix: Proofs of Theorem 1 and Lemma 4 Proof of Theorem 1.

Consider the c-evolution at a generic time t  N to commands g(?t) 2 R, 26

g(t)  g 2 R, 8t 2 Z+, c(t) =

?1 X k=?1

^ 1g + W~ 1 g ? Ht?k g(k) + W

1 X k=t+1

Ht g

Since R is bounded, all references r 2 R satisfy inequalities of the form jrj j  rj , 8j = 1; : : : ; p, where rj denotes the j -th component of r. Therefore, for t  N ? 1, ~ 1 g]i j  2  X M ij rj 1 +W 1 ?  j =1 t+1 p

jci (t) ? [W^ By setting

P

log[ 2 (1 ? ) pj=1 M ij rj ] ti = ?1 log 

(32)

t = max ti i=1;:::;q

(33)

jci (t) ? c^ig ? [W~ 1g]i j  ; 8t  t; 8i = 1; : : : ; q

(34)

and

one has

Since g 2 G , it follows that c(t) 2 C , 8t > t . This proves the \)" part. The \(" part is obvious. 2

Proof of Lemma 4.

Without loss of generality, assume  = 0. De ne g~(t) , g(t) ? g0 and consider the predicted evolution of vector c obtained by supplying system  with the constant command g(t) = g , 8t 2 Z+ c(t) =

?1 X k=?1

Ht?k [g0 + g~(k)] +

Xt k=0

Ht?k g

27

= (W^ 1 + W~ 1 )[g0 + (g1 ? g0 )] + ?1 X

k=?1

Ht?k [~g(k) ? (g1 ? g0 )]

(35)

By Lemma 3, g 2 G . In order to prove that g is admissible, we must show that c(t) 2 C , 8 2 S , and 8t  0, or, equivalently, by (35), P?k=1?1 Ht?k [~g(k)?

(g1 ? g0)]  , 8t  0. Then, ?1 X k=?1

Ht?k [~g (k) ? (g1 ? g0 )] =

1 X

k=maxfN;t+1g

Hk g~(k) +

1 X k=t+1

NX ?1 k=t+1

Hk g~(k) +

Hk (g1 ? g0 )

where the rst sum in the second term is equal to zero for t  N ? 1. By letting

9 8 8 > >> ^ i >> > >> jWN ?1 ? W^ tj+ > >> > = >> max < ij ij ? W W j ; max fj >> N ?1 t >; > > > < i = 1; : : : ; q >> ij > j ij > 1 , t = 0; : : : ; N ? 2 : jW N ?1 ? W t jg ; > >> >> if 0  t  N ? 2 > >: 0; if t  N ? 1 N=2 ij E i=1;:::;q 1 ? 

j2 , max

it follows that c(t)  Bc ?  +

p 2 X X h=1 j =1

28

jh[ + (g1j ? g0j ) ]  Bc

for

8 >>    p >> 2 X X >< 2 jh h=1 j =1 >>  >>  X p 2 X >: 2 jh (g1j ? g0j ) h=1 j =1

(36)

2

29

List of Figures 1

Step-response interval ranges

6

2

Impulse-response interval ranges

6

3

Reference set

8

4

Servomechanism model

18

5

Unconstrained linear response

21

6

Uncertainty ranges

22

7

Response for JL = 25JM

23

8

Response for di erent uncertainty ranges

24

9

Fast controller + PRF, and linear controller

25

List of Tables 1

Model parameters

19

30