Model Predictive Control for Nonlinear Sampled-data Systems
L.Gr¨une1, D.Neˇsi´c2, and J.Pannek3
1 Mathematical Institute, University of Bayreuthlars.gruene@uni-bayreuth.de
2 EEE Department, University of Melbourne, Australiad.nesic@ee.mu.oz.au
3 Mathematical Institute, University of Bayreuth juergen.pannek@uni-bayreuth.de
Summary. The topic of this paper is a new model predictive control (MPC) ap- proach for the sampled–data implementation of continuous–time stabilizing feedback laws. The given continuous–time feedback controller is used to generate a reference trajectory which we track numerically using a sampled-data controller via an MPC strategy. Here our goal is to minimize the mismatch between the reference solution and the trajectory under control. We summarize the necessary theoretical results, discuss several aspects of the numerical implemenation and illustrate the algorithm by an example.
Key words: nonlinear control, sampled-data, model predictive control, reference solution
1 Introduction
Instead of designing a static state feedback with sampling and zero order hold by designing a continuous–time controller which is stabilizing an equi- librium and discretizing this controller ignoring sampling errors which leads to drawbacks in stability, see [5, 8], our approach is to use a continuous–time feedback and to anticipate and minimize the sampling errors by model pre- dictive control (MPC) with the goal of allowing for large sampling periods without loosing performance and stability of the sampled–data closed loop.
Therefore we consider two systems, the first to be controlled by the given continuous–time feedback which will give us a reference trajectory, and a sec- ond one which we are going to control using piecewise constant functions to construct an optimal control problem by introducing a cost functional to measure and minimize the mismatch between both solutions within a time interval.
In order to calculate a feedback instead of a time dependent control function and to avoid the difficulties of solving a Hamilton-Jacobi-Bellman equation for an infinite horizon problem we reduce the infinite time interval to a finite one by introducing a positive semidefinite function as cost–to–go. To re–gain the infinite control sequence we make use of a receding horizon technique.
For this approach we will show stability and (sub-)optimality of the solution under certain standard assumptions.
We will also show how to implement an algorithm to solve this process of iteratively generating and solving optimal control problems. The latter one is done using a direct approach and full discretization that will give us one optimization problem per optimal control problem which can be solved using an SQP method.
Therefore in Section 2 the problem, the necessary assumptions and our control scheme will be presented. In Section 3 we review the theoretical background results about stability and inverse optimality from [13]. Having done this the numerical implementation will be presented and discussed in Section 4 and its performance will be demonstrated by solving an example in Section 5. Finally conclusions will be given in Section 6.
2 Problem formulation
The set of real numbers is denoted asR. A functionγ:R≥0→R≥0 is called classG if it is continuous, zero at zero and non-decreasing. It is of classK if it is continuous, zero at zero and strictly increasing. It is of classK∞ if it is also unbounded. It is of classLif it is strictly positive and it is decreasing to zero as its argument tends to infinity. A functionβ:R≥0×R≥0→R≥0 is of classKLif for every fixedt≥0 the functionβ(·, t) is of classK and for each fixeds >0 the functionβ(s,·) is of classL. Given vectorsξ, x∈Rn we often use the notation (ξ, x) := (ξT, xT)T.
We consider a nonlinear feedback controlled plant model
˙
x(t) =f(x(t), u(x(t))) (1)
with vector fieldf :Rn×U→Rn and statex(t)∈Rn, whereu:Rn→U⊂ Rm denotes a known continuous–time static state feedback which (globally) asymptotically stabilizes the system. We want to implement the closed loop system using a digital computer with sampling and zero order hold at the sampling time instants tk =k·T, k∈N, T ∈R>0. Then for a feedback law uT(x) the sampled-data closed loop system becomes
˙
x(t) =f(x(t), uT(x(tk))), t∈[tk, tk+1). (2) Our goal is now to designuT(x) such that the corresponding sampled–data solution of (2) reproduces the continuous–time solutionx(t) of (1) as close as possible.
The solution of the system (1) at time t emanating from the initial state x(0) = x0 will be denoted by x(t, x0). Also we will assumef(x, u(x)) to be locally Lipschitz in x, hence a unique solution of the continuous–time closed loop system to exist for any x(0) = x0 in a given compact set Γ ⊂ Rn containing the origin.
Remark 2.1 The simplest approach to this problem is the emulation design in which one simply sets uT(x) := u(x). This method can be used for this purpose but one can only prove practical stability of the sampled–data closed loop system if the sampling timeT is sufficiently small, see [8].
In order to determine the desired sampled–data feedback uT we first search for a piecewise constant control functionv whose corresponding solution ap- proximates the solution of the continuous–time closed loop system. Therefore the mismatch between the solutions of
˙
x(t) =f(x(t), u(x(t))), x(t0) =x0 (3) ξ(t) =˙ f(ξ(t), v[0,∞]), ξ(t0) =ξ0 (4) can be measured. Hereξ(t, ξ0) denotes the solution of the system under control and v[0,∞] is a piecewise constant function with discontinuities only at the sampling instants tk :=k·T, k∈N. In order to measure and minimize the difference between both trajectories a cost functional of the form
J(ξ(t), x(t), v[0,∞)) :=
∞
X
j=0 T
Z
0
l(ξ(t)−x(t), vj)dt (5) is needed wherel:Rn×U→R≥0. This results in an optimal control problem with infinite horizon which involves solving a Hamilton-Jacobi-Bellman type equation. In the linear case solutions to differentH2andH∞control designs are known but the nonlinear case is typically too hard to be solved.
In order to avoid this computational burden we consider a reduced problem in a first step by limiting the horizon to a finite length. This will give us a suboptimal MPC controller whose numerical computation is manageable.
SinceT is fixed due to the problem formulation the length of the horizonH can be given by M ∈ N viaH = M ·T. Hence the cost functional can be written as
JM(ξ(t), x(t), v[0,M−1]) :=
M−1
X
j=0 T
Z
0
l(ξ(t)−x(t), vj)dt +F(ξ(tM), x(tM)) (6) using the function F to measure the cost-to-go
∞
P
j=M T
R
0
l(ξ(t)−x(t), vj)dt.
Remark 2.2 It is not necessary for F to be a control-Lyapunov-function of (3), (4) to prove semiglobal practical stability of the closed loop system.
Moreover terminal costs of the form F(ξ(tM), x(tM)) instead of F(ξ(tM)− x(tM)) are considered since the infinite horizon value function V∞(ξ, x) :=
v[0,∞)inf
J(ξ, x, v[0,∞))does not have in general the formV∞(ξ−x).
Using this approach an optimal control problem with finite horizon has to be solved which will return a finite control sequence ˆu[0,M−1]. In order to deter- mine the sampled–data feedback lawuT an infinite sequence of optimal control problems can be generated and solved using a receding horizon approach. To this end in a second step only the first control
u=uM(ξ, x) := ˆu0(ξ, x) (7) is implemented and the horizon is shifted forward in time byT. Hence a new optimal control problem is given and the process can be iterated. Accord- ing to this procedure the receding horizon control law uT = uM is a static state feedback for the coupled system that is implemented in a sampled-data fashion. Then the overall closed loop system is given by
ξ(t) =˙ f(ξ(t), uM(ξ(tk), x(tk))), ξ(0) =ξ0, t∈[tk, tk+1), (8)
˙
x(t) =f(x(t), u(x(t))), x(0) =x0. (9)
Remark 2.3 We like to emphasize that it is not only our goal to obtain asymptotical stability of (8), (9) which implies tracking since we have that
|(ξ(t), x(t))| ≤β(|(ξ0, x0)|, t) ∀t≥0, (10) but also that we achieve this in an appropriate sub–optimal manner.
3 Stability and Inverse Optimality
Since most of the time one can only work with approximated discrete-time models consistency with the exact discrete-time model as described in [11]
is needed. Under the consistency condition given by Definition 1 in [10] and suitable mild additional assumptions one can conclude that asymptotic stabil- ity of the approximate model carries over to semiglobal practical asymptotic stability for the exact model, see [11] for a general framework and [3] for cor- responding results for MPC algorithms. This justifies the use of numerical approximations, cf. also Remark 4.1, below. To conclude semiglobal asymp- totical stability of the closed loop system using the proposed MPC controller we present the following theorem, which relies on Theorem 1 in [1].
Theorem 3.1 (Stability)
Suppose the following conditions hold:
1. l andF are continuous;
2. Uis bounded;
3a. The continuous–time system (1)is globally asymptotically stable;
3b. There exists a constantr0>0 and a function γ∈ K∞ with l(y, u)≥max
max
|x|≤2|y||f(x, u)|, γ(|y|)
, ∀|y| ≥r0; 3c.f(·,·)andu(·)are locally Lipschitz in their arguments;
4. The value function is such that for some α¯ ∈ K∞ we have that Vi(ξ, x)≤ α(|(ξ, x)|)¯ for all i≥0 and all(ξ, x)∈R2n.
Then there exists a functionβ∈ KLsuch that for each pair of strictly positive real numbers(∆, δ)there exists a constantM1∗∈Z≥1such that for all(ξ, x)∈ B∆ andM ≥M1∗ the solutions of the continuous–time system(8),(9)satisfy
|(ξ(t), x(t))| ≤max{β(|(ξ0, x0)|, t), δ} ∀t≥0. (11) Proof. Make use of the underlying discrete-time system via Theorem 1 in [1]
and Theorem 1 in [10], see [13] for details.
Therefore one can apply the calculated MPC control coming out of an ap- proximated model in reality without loss of stability.
Remark 3.2 IfF is a control Lyapunov function for the exact discrete-time model of the uncontrolled sampled-data system
ξ+ =G(ξ, u) :=ξ(T, ξ, u), ξ(0) =ξ0, (12)
x+ =H(x) :=x(T, x), x(0) =x0 (13)
then it follows from [1] that the theoretical bound M1∗ for the necessary pre- diction horizon decreases which was confirmed in our numerical simulations.
Remark 3.3 Explicit bounds to guarantee the stability properties of the un- derlying discrete-time system can be found in [1].
In order to show inverse optimality of our approach we suppose thatFis such that there exists a closed setXf ⊂R2n and a control lawu=uf(ξ, x) with 1. uf(ξ, x)∈U ∀(ξ, x)∈Xf
2. If (ξ, x) is∈Xf then also (G(ξ, uf(ξ, x)), H(X)) is∈Xf. 3. For all (ξ, x)∈Xf we have that
F(G(ξ, uf(ξ, x)), H(x))−F(ξ, x)≤ −
T
Z
0
l(ξ(s, ξ, uf)−x(s, x), uf)ds.
Theorem 3.4 (Inverse (Sub-)Optimality)
Consider the discrete-time model (12), (13) and suppose that the previous assumptions are valid. Then there exists a set XM ⊂ R2n and a function Q:Rn×Rn×U→Rwith
Q(ξ, x, uM)≥Q(ξ, x, uM) :=
T
Z
0
l(ξ(s, ξ, uf)−x(s, x), uf)ds (14)
∀(ξ, x)∈XM,u∈Usuch that for all(ξ, x)∈XM we have that the controller (7) minimizes the cost functional
J(ξ, x, u[0,∞)) :=
∞
X
i=0
Q(ξi, xi, ui). (15) Proof. The principle of optimality and the stated assumptions are utilised to showQ(ξ, x, uM)≥Q(ξ, x, uM), see [13] for details.
4 Numerical Solution
For the solution of the optimal control problem we use a direct approach and therefore replace the problem to minimize (6) with dynamics (3), (4) by numerical approximations ˜ξ(t, ξ0, u) ofξ(t, ξ0, u) and ˜x(t, x0, u) ofx(t, x0), respectively. For this approach convergence has been investigated in [7] and under suitable conditions one can guarantee that the order of convergence is O(T).
From this formulation one obtains an optimization problem by introducing the variable z = (ξ0, . . . , ξM, x0, . . . , xM, u0, . . . , uM) and rewriting the ap- proximated optimal control problem as
MinimizeF(z) :=
M−1
X
j=0 T
Z
0
l( ˜ξ(s, ξj, vj)−x(s, x˜ j), vj)ds+F(ξM, xM)
s.t.G(z) :=
[−ξj+1+ ˜ξ(h, ξj, vj)]j=0,...,M−1
[−xj+1+ ˜x(h, xj)]j=0,...,M−1
(ξ0, x0)−(ξ0, x0)
= 0
with the constraints coming along with the approximation. This is a well known problem that can be solved using the KKT conditions by SQP methods if the cost functional and the constraints are sufficiently often differentiable in a sufficiently large neighborhoodN(z∗) of the local minimaz∗. These methods are known to be stable and efficient even for large scale systems.
The used algorithm computes a sequence (z[k]) via z[k+1] = z[k] +α[k]p[k]. Within this iteration the search directionp[k] is calculated by generating and solving quadratic subproblems of the form
min
p∈RNz
∇zF(z[k])p+1
2pTB[k]p s.t. G(z[k]) +∇zG(z[k])p= 0.
The algorithm computes the derivatives by forward difference schemes if they are not given by the user and the matrix B[k] is an approximation of the Hesse matrix where a BFGS-Rank 2 update is implemented so that the Hesse matrix has to be calculated only once. Therefore the usual quadratic order of convergence of the Newton method is reduced but superlinear convergence can still be shown. The step sizeα[k] is obtained by minimizing a merit function L(z, η, ρ) =˜ L(z, η) + 12
Nz
P
j=1
ρjG2j(z) such that the natural step size α[k] = 1 of the Newton method is reduced but one can expect it to be close to 1 in a small neighborhood ofz∗.
Remark 4.1 Since our aim is to allow for large sampling periodsT an adap- tive step size control algorithm such as DoPri5, see [4], is necessary within each interval [kT,(k+ 1)T),k∈N, in order to avoid errors in the state tra- jectories and the cost functional which therefore has to be transformed. Note that the local accuracy guaranteed by the step size control here plays the role of the accuracy parameterδ in the consistency Definition 1 in [10].
Remark 4.2 The case of an integration step sizehthat is different from the sampling periodT has been analysed theoretically for MPC schemes in [3]. An important aspect of this analysis is that h, or — more generally — the nu- merical accuracy parameter, can be assigned arbitrarily and independently of T (where of course one has to ensure that the sampling instants are included in the set of gridpoints used for integration in order to match the discontinu- ities of the control function). It should be noted that our algorithm fulfils this requirement. In fact, when we remove the x–subsystem (3) from our scheme and use a local Lyapunov function as a terminal cost we obtain exactly the direct MPC algorithm discussed theoretically in [3].
Compared to this standard MPC approach the main difference of our scheme lies in the fact that we can directly enforce a desired transient behavior in- duced by the continuous–time feedback, while in standard MPC schemes the transient behavior can only be influenced indirectly through the choice of the cost functional. Clearly, the design of a suitable continuous–time controller requires a considerable amount of a priori work, but this may be rewarded by a better performance of the resulting sampled–data closed loop.
Remark 4.3 An important problem is the choice of a good initial guess v[0,M−1] for the optimization, keeping in mind that we deal with a nonlin- ear optimization problem. Even though suboptimal solutions to this problem may be sufficient to ensure stability, see [6], here we also aim at good per- formance. Convergence to the global optimum, however, can only be expected
when the initial solution is already close to it. When passing from tk totk+1 the shifted optimal control sequence from the previous step typically yields such a good initial guess, which is confirmed by our numerical experience that the computational costs for the optimization are decreasing monotonically during the iteration process.
A more severe problem is the choice of the initial guess att0when no previous optimal control is known. In this case, in our approach the known continuous–
time feedback can be exploited for this purpose when the emulated feedback from Remark 2.1 yields solutions which do not deviate too far from the continuous–
time reference. However, this method fails when the emulated feedback leads to unstable solutions and the time horizon H = M ·T is rather large. Such situations can sometimes be handled by reducing the length of the horizon H =M·T but proceeding this way one has to keep in mind that there exists a lower bound forH from the stability proof. Also, simulations have shown that while on one hand computational costs grow with the length of the horizon, on the other hand better performance can be achieved using longer horizons.
Therefore, at the moment it is up to simulations to find a good set of param- eters and a good initial guess ofv[0,M−1].
A promising alternative approach and topic of future research is whether some of the methods developed in [9, 12] can be used in order to construct the initial guess, an approach that would lead to a predictor–corrector type algorithm in which the MPC strategy plays the role of the corrector.
5 Example
Here we present a model of a synchronous generator taken from [2]
˙
x1=x2, x˙2=−b1x3sinx1−b2x2+P, x˙3=b3cosx1−b4x3+E+u. (16)
0.2 0.7
1.2
−2
−1 0 1 2 1 1.5 2
Emulation MPC
0.2 0.7
1.2
−2
−1 0 1 2 1 1.2 1.4 1.6 1.8 2
Emulation MPC
0.2 0.7
1.2
−2
−1 0 1 2 1 1.5 2
Reference
Fig. 1.Improvement by MPC control over emulation fora1= 0.45,T = 0.1 (left), T = 0.5 (middle) and reference solution using continuous–time feedback (right) We use the parameter b1 = 34.29, b2 = 0.0, b3 = 0.149, b4 = 0.3341, P = 28.22 and E = 0.2405, as well as the continuous–time feedback law u(x) = a1((x1−x∗1)b4+x2) with feedback gaina1>0, whose purpose is to enlarge the domain of attraction of the locally stable equilibrium x∗ ≈ (1.12,0.0,0.914) (note that this equilibrium is locally asymptotically stable also foru≡0). As initial value we used the vector x0 = (0.5,0.0,2.0) and generated results for T = 0.1,T = 0.5 anda1= 0.45.
One can see that the fast dynamics of the problem require small sampling periods to maintain stability using the emulated feedback law. The MPC
control on the other hand not only stabilizes the equilibrium even for rather largeT but also keeps the sampled–data solution close to the reference.
6 Conclusion
We proposed an unconstrained model predictive algorithm for the sampled–
data implementation of continuous–time stabilizing feedback laws. Stability and inverse optimality results were briefly revisited and numerical issues were discussed. Compared to direct MPC approaches without using continuous–
time feedbacks, advantages of our method are that the sampled–data solutions inherit the performance properties of the continuous–time controller and that the knowledge of the continuous–time controller helps to reduce the com- putational cost of the numerical optimization. Future research will include a systematic study about how this knowledge can be used in a numerically efficient way and an extension of our approach to dynamic continuous–time controllers.
References
1. Grimm G, Messina M J, Teel A R, Tuna S (2005) Model predictive control:
for want of a local control Lyapunov function, all is not lost. In: IEEE Trans.
Automat. Contr. 50: 546–558
2. Gr¨une L (2001) Subdivision Techniques for the Computation of Domains of Attractions and Reachable Sets. In: Proceedings of NOLCOS 2001: 762–767 3. Gyurkovics E, Elaiw A M (2004) Stabilization of sampled–data nonlinear sys-
tems by receding horizon control via discrete–time approximations. In: Auto- matica 40: 2017–2028
4. Hairer E, Nørsett S P, Wanner G (1993) Solving Ordinary Differential Equa- tions I. Nonstiff Problems. 2nd Ed., Springer–Verlag
5. Jadbabaie A, Hauser J, Yu J (2001) Unconstrained receding horizon control of nonlinear systems. In: IEEE Trans. Automat. Contr. 46: 776–783
6. Jadbabaie A, Hauser J (2005) On the stability of receding horizon control with a general terminal cost. In: IEEE Trans. Automat. Contr. 50: 674–678 7. Malanowski K, B¨uskens Ch, Maurer H (1997) Convergence of approximations
to nonlinear optimal control problems. In Fiacco, A. (ed.) Mathematical Pro- gramming with Data Perturbations. Marcel Dekker
8. Laila D S, Neˇsi´c D, Teel A R (2002) Open and closed loop dissipation inequali- ties under sampling and controller emulation. In: Europ. J. Control 18: 109–125 9. Monaco S, Normand–Cyrot D (2001) Issues on nonlinear digital control. In:
Europ. J. Control 7: 160–178
10. Neˇsi´c D, Teel A R, Sontag E D (1999) Formulas relating KL stability estimates of discrete–time and sampled–data nonlinear systems. In: Syst. Contr. Lett. 38:
49–60
11. Neˇsi´c D, Teel A R (2004) A framework for stabilization of nonlinear sampled–
data systems based on their approximate discrete–time models. In: IEEE Trans.
Automat. Contr. 49: 1103-1122
12. Neˇsi´c D, Gr¨une L (2005) Lyapunov based continuous–time nonlinear controller redesign for sampled–data implementation. In: Automatica 41: 1143–1156 13. Neˇsi´c D, Gr¨une L (to appear) A receding horizon control approach to sampled–
data implementation of continuous–time controllers. In: Syst. Contr. Lett.
