Model predictive sampled-data redesign for nonlinear systems

We propose a model predictive control (MPC) strategy for sampled-data implementation (with the zero order hold assumption) of continuous-time controllers for general nonlinear systems. We assume that a continuous-time controller has been designed so that the continuous-time closed-loop satisfies all performance requirements. Then, we use this control law indirectly to compute numerically a sampled-data controller via an MPC strategy that minimizes the mismatch between the solutions of the sampled-data model and the continuous-time closed-loop model. We present conditions under which stability and sub-optimality of the closed loop can be proved.


I. INTRODUCTION
Nowadays, modern controllers are typically implemented digitally and this fact strongly motivates investigation of sampled-data systems that consist of a continuous-time plant controlled by a discrete-time (digital) controller.While tools for analysis and design of linear sampled-data systems are well developed, similar results for nonlinear systems still need development to be as useful as their linear counterparts.A possible approach for sampled-data controller design is to first design a continuous-time controller for the continuoustime plant ignoring sampling and then discretize the obtained controller for digital implementation [3], [5], [12].The classical discretization methods, such as the Euler, Tustin or matched pole-zero discretization are attractive for their simplicity but they may not perform well in practice since the required sampling rate may exceed the hardware limitations even for linear systems [1], [10].This has lead to a range of advanced controller discretization techniques based on optimization ideas that compute "the best discretization" of the continuous-time controller in some sense.A nice account of these optimization based approaches for linear systems has been given in the Bode Lecture by Anderson in [1] and later in the book [3].
We are not aware of a similar optimization based approach for discretization of continuous-time controllers for nonlinear systems.A possible reason for this may be the inherent computational complexity of nonlinear optimal control problems that inevitably require solutions to Hamilton-Jacobi type equations.However, while nonlinear optimal controllers are often impossible to compute in practice due to the computational burden associated with solving the Hamilton-Jacobi equations, different suboptimal solutions are much more tractable.For example, the model predictive (or receding horizon) control has a manageable computational complexity for relatively large nonlinear problems [7], [9], [14].
It is the purpose of this paper to describe and numerically illustrate a novel model predictive control scheme that can be used in implementing digitally continuous-time controllers that have been already designed.The cost function that we consider penalizes the difference of solutions of the continuous-time closed-loop system and the sampled-data solutions.In this sense, the control scheme that we consider can be regarded as a nonlinear and sub-optimal version of linear results presented in [1], [3].For simplicity, we only consider digital implementation of static state feedback controllers and we present results only for unconstrained MPC.Stability of our control scheme follows from recent results in [7] and is analyzed in detail in [15].Under appropriate conditions our MPC scheme is also inverse optimal (or sub-optimal) in some sense, which is similar to standard MPC results, see [14].
The paper is organized as follows.In Section II we present preliminaries, pose the problem we consider and present our control scheme.Details about its stability and inverse optimality properties are given in Section III.Its numerical implementation and numerical examples are presented in Section IV.The final section V concludes our paper.

II. PRELIMINARIES
The set of real numbers is denoted as strictly positive and it is decreasing to zero as its argument tends to infinity.A function β : R ≥0 × R ≥0 → R ≥0 is of class KL if for every fixed t ≥ 0 the function β(•, t) is of class K and for each fixed s > 0 the function β(s, •) is class L. Given vectors ξ, x ∈ R n we often use the notation (ξ, x) := (ξ T x T ) T .

A. Problem formulation
Consider the plant: where x ∈ R n and u ∈ U ⊂ R m are respectively the state and the control input of the system.Standing assumptions that we will use throughout the paper are as follows: Standing Assumptions: 1) A continuous-time controller u = u(x), such that u(x) ∈ U, ∀x ∈ R n has been designed for the continuous-time plant (1) so that the continuous-time closed-loop system is (globally) asymptotically stable and satisfies all performance requirements.2) The controller is to be implemented using a sampler and zero order hold.In other words, for a given fixed sampling period T > 0 the control signal is constant during sampling intervals, i.e. u(t) = u(t k ) = const., ∀t ∈ [t k , t k+1 ), k ∈ N, where t k := kT .We will always use x(t, x 0 ) to denote the solution of the system (2) at time t emanating from the initial state x(0) = x 0 .We assume in what follows that f (x, u(x)) is locally Lipschitz in x and, hence, for any x(0) = x 0 the continuous time closed loop system (2) has a unique solution.
Before we pose the problem we consider, we first define what we mean by solutions of the sampled-data system.This definition is the same as the definition of S-solution (sampled solution) proposed in [4].Given an initial state ξ(t 0 ) = ξ 0 and a control signal v(t) = v k , t ∈ [t k , t k+1 ), k ∈ N, the solution of the sampled data system on interval [t 0 , t 1 ] is the solution of the continuous time system: Let the solution of this system at time t 1 = T be denoted as ξ(t 1 ).Then, solution of the sampled-data system on the time interval [t 1 , t 2 ] is the solution of the continuous-time system and so on.Denote a sequence of controls The solution of the sampled data system at time t, starting at ξ 0 and under the sequence v [0,M ] is denoted as ξ(t, ξ 0 , v [0,M ] ) or simply ξ(t) when ξ 0 and v [0,M ] are clear from the context.The problem that we consider in this paper is as follows: Find a sampled data controller so that for any given ξ(0) = ξ and x(0) = x, the solution ξ(t) of the sampled-data system reproduces the solution x(t) of the continuous-time system "as close as possible".The optimal solution of the above problem would necessarily involve infinite horizon minimization (optimization) of some measure of the mismatch between sampled-data and continuous-time solutions (ξ(•) and x(•)).While for linear systems this problem is feasible, cf.[1], [3], for nonlinear systems infinite horizon optimization typically leads to computationally intractable problems: Hence we will instead investigate suboptimal model predictive (or receding horizon) controllers that are known to be much more manageable computationally.In particular, given a fixed positive integer M , the controllers we consider involve minimization (in ) of the cost of the form: In particular, at each sampling interval we solve the following unconstrained optimization problem: where the actual plant state ξ is measured at sampling times and x is the nominal reference state from (2), which is determined numerically.Moreover, we implement the controller in a receding horizon fashion where at each sampling interval we apply only the first control in the optimal sequence û[0,M−1] .At the next sampling interval the new control sequence is obtained by solving again the optimization problem with new measured states and only the first control in the sequence is actually applied.Note that the receding horizon control law is a static state feedback u = u M (ξ, x) that is implemented in a sampleddata fashion so that the overall closed loop system can be written as follows: with ξ(0) = ξ and x(0) = x.Remark 2.1: One goal of our analysis will be to show that the system (7) is asymptotically stable since then we have for some β ∈ KL, which implies tracking, i.e. lim t→∞ |ξ(t) − x(t)| = 0.It is obvious from item 1) of our Standing Assumptions that any sampled-data controller that yields lim t→∞ |ξ(t)| = 0 would imply tracking.However, model predictive controllers that we propose are sub-optimal in an appropriate sense (see Theorem 3.5) and, hence, they are not only achieving tracking but also do so in an appropriate sub-optimal manner.Remark 2.2: Note that in the system (7) the control designer can choose to initialize the bottom continuous-time subsystem in a particular manner, since this is just a reference model to be used in computing the controller u M (ξ, x).For instance, we could measure the initial state of the sampleddata system at the initial time ξ(0) = ξ and then let x(0) = ξ(0) = ξ.This makes sense because we would like the sampled-data system to recover as close as possible behavior of the continuous-time system from the same initial condition.Nevertheless, we will present analysis of stability of the system (7) that yields bounds on transients for arbitrary initializations.
Remark 2.3: A brute force approach to solving the above problem is to simply implement the emulated controller: and then sample as fast as possible (reduce T ).This approach was shown in [12] to recover the performance of the continuous-time system in an appropriate semi-global practical sense (T is the parameter that we need to reduce sufficiently).However, due to hardware limitations on the minimum achievable T this approach is often not feasible.
In this paper, we use the designed continuous-time control law indirectly as a part of a receding horizon strategy.Remark 2.4: Note that we consider general terminal costs of the form F (ξ(t M ), x(t M )) instead of the costs of the form F (ξ(t M ) − x(t M )) that may appear to be more natural in this context.However, if we think of the terminal cost as the approximation of the infinite horizon value function, then it is obvious that the form of F (•, •) that we use is more appropriate since the infinite horizon value function would not in general have the form V ∞ (ξ − x).
Remark 2.5: The problem we consider can be viewed as a special case of a more general problem where we assume that the reference model which we would like our sampleddata closed-loop to track is not necessarily of the form (2) but perhaps of the more general form where we assume that this reference model is stable and satisfies all performance requirements.

B. Stability properties of discrete-time MPC schemes
Our stability results will be heavily based on recently proved stability results of discrete-time MPC schemes in [7].A unique feature of these results is that the terminal and stage costs do not have be positive definite functions of the state, which is the case for our MPC scheme.Moreover, the terminal cost does not have to be a local control Lyapunov function for the system to show stability.In this section, we summarize results from [7] and in the next section we show how they can be used to analyze stability of the closed loop system (7).To this end, we introduce an auxiliary discretetime problem since [7] deal only with discrete-time systems.Note that given any fixed sampling period T > 0 the (exact) discrete-time model of the uncontrolled sampled-data system (3), (2), when it is well defined, has the form: where G(ξ, u) := ξ(T, ξ, u) and H(x) := x(T, x).We denote the solutions of the discrete-time model (9) as ξ(k, ξ, v [0,k−1] ) and x(k, x) or simply by ξ k and x k when the initial states and the control sequence are clear from the context.Moreover, by introducing the following: we can rewrite the cost (4) as follows: where Q is defined in (10) and ξ i , x i are solutions of the discrete-time system (9).Finally, optimization problem (5) can be rewritten as follows: and the discrete-time model of the closed-loop sampled-data system (7) can be written as follows: where u M (•, •) comes from (6).Remark 2.6: It is a standard result in the literature to show (under weak assumptions) that stability of the discrete-time model (13) implies stability of the sampled-data system (7) (see, for instance [18]).
Remark 2.7: Typically the (exact) discrete-time model of the system is not available and hence the controller design needs to be based on an approximate discrete-time model.We do not investigate these issues in this paper and refer to the rigorous framework for controller design based on approximate discrete-time models developed in [17], [16].
We next adapt appropriate definitions and results from [7] to be applicable for stability analysis of (13).Some of these definitions can be further relaxed (see [7]) but the version we present suffices for our purposes.
Definition 2.8: Consider the system (9) and a function Q = Q(ξ, x, u).The system (9) is said to be detectable from Q with respect to (α W , α W , γ W ) if α W , γ W ∈ K ∞ and α W ∈ G and there exists a continuous function W : R 2n → R ≥0 such that for all (ξ, x) ∈ R 2n and all u ∈ U: ) Definition 2.9: We will say that the terminal cost F is a control Lyapunov function for the system if it can be decomposed as , where M is the optimization horizon, the function Γ : Z ≥1 → R ≥1 is nondecreasing and unbounded and, moreover, there exist functions α and for every (ξ, x) ∈ R 2n there exists u ∈ U such that The following results were proved in [7]: Theorem 2.10: Suppose that: (i) Q and F are continuous; (ii) U is bounded; (iii) The system (9) is detectable from Q(ξ, x, u) for some functions (α W , γ W , α W ); (iv) The value function is such that for some α ∈ K ∞ we have that V i (ξ, x) ≤ α(|(ξ, x)|) for all i ≥ 0 and all (ξ, x) ∈ R 2n .Then, for each M ≥ 2 there exist α Y , α Y , α Y ∈ K ∞ , β Y ∈ KL and a continuous function Y M : R 2n → R ≥0 such that for all (ξ, x) ∈ R 2n we have: Moreover, if the following condition also holds: (v) The terminal cost F is a control Lyapunov function for the system (9), then there exists a continuous function Y M and β Y ∈ KL such that for all (ξ, x) ∈ R 2n (16) holds and Remark 2.11: Explicit formulas for computing all bounding functions in Theorem 2.10 are given in [7].
A direct consequence of Theorem 2.10 is Proposition 2.12: Suppose that items (i), (ii), (iii) and (iv) of Theorem 2.10 hold.Then for each pair of strictly positive real numbers (∆, δ) there exists M * 1 ∈ Z ≥1 such that for all (ξ, x) ∈ B ∆ the solutions of the system (13) satisfy: Suppose, moreover, that condition (v) of Theorem 2.10 holds.Then, for each pair of strictly positive real numbers (∆, δ) there exists M * 2 ∈ Z ≥1 such that for all (ξ, x) ∈ B ∆ the solutions of the system (13) satisfy (18).
Remark 2.13: Significance of the condition (v) in the above theorem is that the prediction horizon M * 2 obtained from the proof is typically smaller than M * 1 for the same given (∆, δ).
Remark 2.15: With minor changes one can state regional stability results instead of the global results presented here.

III. STABILITY AND OPTIMALITY
Under reasonable assumptions we can conclude stability of the sampled-data system (7).The following theorem is the main result in this context.
Theorem 3.1: Suppose that the following conditions hold: (i) and F are continuous; (ii) U is bounded; (iii-a) The continuous-time system (2) is globally asymptotically stable; (iii-b) There exists r 0 > 0 and γ ∈ K ∞ with1 l(r, u) ≥ max max The proof consists of showing that conditions of Theorem 3.1 imply that all conditions of Proposition 2.12 hold for the underlying discrete-time system, details can be found in [15].
We note that Remarks 2.13 and 2.14 hold accordingly for this theorem.
Remark 3.2: Note that although we are considering an unconstrained optimization problem, the problem may not be feasible for all (ξ, x) ∈ R 2n and u ∈ U. First, controllability of the sampled-data (and, hence, discrete-time) system can be lost due to sampling.Moreover, due to possible finite escape times, the solutions of the sampled-data system may not be defined for all initial states and inputs.However, it is standard to show that both of these issues disappear on any given compact set if we sample fast enough.
Remark 3.3: Note that our problem can be interpreted as a tracking problem where we want trajectories of the sampleddata system to track the trajectories of the continuous-time model.This set up is standard in the output tracking literature where the continuous-time model is termed an exogenous model (e.g.see [8]).We note, however, that our controller does not have the typical internal model structure because our exogeneous system is not Poisson stable.
It is interesting to note that similar to standard results in model predictive control literature (see [14]), under extra assumptions our control law is inverse optimal.
Assumption 3.4: Suppose that F is such that there exists X f ⊂ R 2n and a control law u = u f (ξ, x) such that the following conditions hold: . For details on how to construct F and meeting Assumption 3.4 we refer to [15].
The following theorem is our main result on inverse optimality.
Theorem 3.5: Consider the discrete-time plant model (9).Suppose that Assumption 3.4 holds.Then, there exists a set Proof: By checking the principle of optimality one proves the optimality for Q(ξ, x, u) [15] for details.

IV. NUMERICAL IMPLEMENTATION AND EXAMPLE
In order to calculate a numerical solution of the optimal control problem a direct discretization method is used.
Following [13], the optimal control problem can be replaced approximatively using numerical approximations ξ(t, ξ, u) of ξ(t, ξ, u) and x(t, x, u) of x(t, x), respectively, to compute arg inf v [0,M −1] J 0 (ξ, x, v [0,M −1] ), where Using the Euler approximation one can guarantee that the order of convergence of this approximation is O(T ), see [13].However, in our MPC scheme the purpose is to use a large sampling rate T which would cause large errors using Euler.To avoid this conflict we generate ξ and x by the fifth order Runge-Kutta-England method with adaptive step size control on each interval [jT, (j + 1)T ], j = 0, . . ., M − 1.The same scheme is used for the numerical approximation of the integral appearing in the definition of J 0 .
The optimization problem is now obtained by introducing the variable z := (ξ 0 , . . ., ξ M , x 0 , . . ., x M , u 0 , . . ., u M ) ∈ R Nz , N z = (2n + m) • (M + 1) to rewrite the approximated optimal control problem as min z∈R Nz F (z) with and z satisfying the constraints This is a well known problem which can be solved using the KKT conditions by SQP methods.These are known to be stable and efficient even for large scale systems and require the functions F and G to be sufficiently often differentiable in a sufficiently large neighborhood N (z * ) of the local minima z * .The algorithm was first presented in [6] and computes a sequence (z [k] ) via z Within this iteration the search direction p [k] is calculated by generating and solving quadratic subproblems where B [k] is an approximation of the Hesse matrix using a BFGS update 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, see [19].The step size α [k] is obtained by minimizing a merit function 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 of z * .We illustrate the performance of our algorithm by a third order Galerkin approximation of the Moore-Greitzer model A stabilizing continuous time feedback has been obtained in [11] and is given by Using the parameters (c 1 , c 2 , σ) = (1, 2, 2), initial value (φ, ψ, R) = (6, 25, 1) and sampling rate T = 0.05, Figure 1 shows that using emulation of u according to Remark 2.3 stabilizes the system but that the model predictive control algorithm improves on the emulation design since the trajectories of the system stay closer to solutions of the system with continuous feedback.In this case the length of the horizon was chosen to be 10T .Changing the above parameters to (c 1 , c 2 , σ) = (1, 50, 2), it can be seen from Figure 2 that emulation of u does no longer stabilize the system.However, also in this case the implemented MPC scheme is able to generate a control sequence that stabilizes the system and keeps the sampleddata solution close to the continuous time one.Note that we achieved these numerical results without using a Lyapunov function for the terminal cost.