Pseudospectral methods for optimal motion planning of differentially flat systems

The article presents some preliminary results on combining two new ideas from nonlinear control theory and dynamic optimization. We show that the computational framework facilitated by pseudospectral methods applies quite naturally and easily to Fliess' implicit state variable representation of dynamical systems. The optimal motion planning problem for differentially flat systems is equivalent to a classic Bolza problem of the calculus of variations. We exploit the notion that derivatives of flat outputs given in terms of Lagrange polynomials at Legendre-Gauss-Lobatto points can be quickly computed using pseudospectral differentiation matrices. Additionally, the Legendre pseudospectral method approximates integrals by Gauss-type quadrature rules. The application of this method to the two-dimensional crane model reveals how differential flatness may be readily exploited.


B. Additive Cascades
Consider the hybrid system _ x1 = g(x1 x2)x1 _ x 2 = 0 2x 2 _ z = 0z x + 2 = x2 + z with reset times t k = kh where k ranges over the nonnegative integers. This system is the additive cascade of GES subsystems since the z-subsystem is GES and the (x 1 ; x 2 ) subsystem is GES when z = 0.
The following result follows from Proposition 6.

V. CONCLUSION
We have extended the results in [7] by constructing an example where the disturbance can be a simple decaying exponential, and by making explicit that the system can be GES with zero input and can have linear sector growth. This enables making observations about additive cascades of globally exponentially stable systems, and about the gradients of Lyapunov functions for globally exponentially stable systems with linear sector growth. We have also provided similar examples for discrete-time and hybrid systems.
Although checking for flatness is an open issue, a growing number of dynamical systems in engineering have been shown to be flat, (see [3] and the references contained therein). For a flat system, the motion planning problem simply reduces to finding a sufficiently smooth output, t 7 ! y(t), that satisfies the boundary conditions in output space.
In principle, finding such smooth functions is not difficult, since the output can be represented in terms of a polynomial with unknown coefficients. These coefficients can then be determined by imposing the condition that the polynomial should satisfy the boundary conditions; however, when differentiating polynomials, it is extremely important to be cognizant of instabilities like the Runge phenomenon [4], [5] associated with interpolating polynomials at equidistant points. Further, in many applications, particularly those arising in astronautics, it is not enough to find feasible trajectories but optimal trajectories that optimize a scalar cost functional given in a Bolza form. For differentially flat systems, the optimal control problem reduces to a classic unconstrained calculus-of-variations problem [6].
In this note, we show that for a differentially flat system, an optimal smooth output function and its derivatives can be easily obtained by pseudospectral methods [4], [5], [7]. Pseudospectral methods are based on approximating the underlying functions by interpolating polynomials which interpolate these functions at some specially chosen Manuscript  nodes. These nodes are the zeros of orthogonal polynomials (or their derivatives) such as Legendre polynomials (Legendre-Gauss points) or Chebyshev polynomials (Chebyshev-Gauss points). Recently, pseudospectral methods have been used very effectively in solving a wide variety of nonlinear optimal control problems as illustrated in [8]- [13].
Here, we show that exploiting differential flatness provides a new way of solving motion planning problems. Our work is similar in spirit to that of Milam et al. [14] and Petit et al. [15] but we will show that our technique is markedly different yet simpler to implement than their B-spline approach. Further, since differentially flat systems are sufficiently smooth, pseudospectral methods provide "exponential convergence rates" for analytic functions and O(N 0m ) for every m for C 1 functions [5] where N is the order of the interpolating polynomial.
This property, known as spectral accuracy, is essentially an outcome of the low Lebesgue constants [4] for Legendre-Gauss and Chebyshev-Gauss node distributions. Spectral accuracy is particularly important for systems where flat outputs cannot be obtained but an output that partially inverts the dynamics can be found and exploited. In this case, the dynamical constraints can be reduced but not eliminated. Potential convergence problems resulting from discretizing the transformed dynamics are handled well by pseudospectral methods.

II. PROBLEM FORMULATIONS
A "classic" smooth optimal control problem can be stated as follows.
Problem C: Determine the trajectory-control pair, [0; f ] 3 7 ! fx x x 2 N ; u u u 2 N g and possibly the clock times 0 and f , that subject to the classic dynamic constraints and end point constraints For simplicity in presentation, we assume all functions to be C 1 -smooth. According to Fliess et al. [1], the dynamical system described by (2) is differentially flat if there exists a variable y 2 N and a function c(1) such that x x x = a(y; _ y; . . . ; y () ) u u u = b(y; _ y; . . . ; y (+1) ) (5) where and are finite positive integers that denote the number of derivatives of the respective variables. The variable y is called a flat or linearizing output. For ease of notation, we let s = +1 and denote the flat output and its derivatives by the composite variable, z 2 (s+1)N z = [y; _ y; . . . ; y (s) ] T (6) so that a(1) : z ! x x x and b(1) : z ! u u u. For a differentially flat system, Problem C can now be replaced by the following.

III. LEGENDRE PSEUDOSPECTRAL METHOD
For the purpose of clarity and brevity, we discuss only the Legendre pseudospectral (PS) method. Let LN(t) be the Legendre polynomial of degree N on the interval [01; 1]. In the Legendre pseudospectral method, the Legendre-Gauss-Lobatto (LGL) points [7], t l , l = 0; . . . N are used. These points are given by t0 = 01, t N = 1; and for 1 l N 0 1, t l are the zeros of _ L N ; the derivative of the Legendre polynomial, L N : For Problem C, the Legendre pseudospectral method offers an approximation for evaluating the integral by Gauss quadratures while the differential constraint is approximated by driving the residuals to zero at the LGL points. In this manner, the Legendre PS method unifies discretization of both the integrals and the derivatives, and in both cases the discretizations are highly accurate. Further details on the approximation method for Problem C are described in [8], [12], and [16]. Here, we focus our attention to Problem DF and the transformations necessary to cast Problem C to this format.
Since the LGL node points lie in the computational interval [01; 1], in the first step of this method, the following affine transformation is used to scale the domain, [0; f ], = (( f 0 0)t + ( f + 0))=2: Next, the vector-valued function, t 7 ! y(t), is written as some Nth where as before the superscript s denotes the sth derivative. It is apparent that we must choose N s + 1. Evaluating the derivatives at t k results in a matrix multiplication of the following form: otherwise.
It can be shown that Di = D i where the superscript denotes matrix powers. Thus, D 2 is obtained by simply squaring D 1 , while D 3 = where B is obtained from (14) and (15). Recall that N is a design parameter and must be chosen such that N s + 1.
In any case, the optimal motion planning problem requires that the integral in (7) be evaluated in terms of the values of the flat outputs and its derivatives at the LGL points. While other polynomial approximations [14] can only use low-order quadrature schemes, in pseudospectral methods, high-order quadrature rules such as the Gauss-Lobatto integration rule can be naturally employed. The integral (7) is approximated by a finite sum which is exact for integrands which are polyno- where w k are the LGL weights [7] w k := 2 N (N + 1) Thus, Problem DF is discretized by the following mathematical programming problem.
IfF andẼ are linear in z k , then for fixed clock times the problem reduces to a linear programming problem for linear conditions in output space. In general, this is a nonlinear programming problem which can be solved using commercial off-the-shelf packages like SNOPT [17]. It is worth noting that in our method, the original state and control variables can be easily recovered by using the differentiation matrices and the functions a(1) and c(1).

IV. EXAMPLE: THE CRANE PROBLEM
A two-dimensional state model of a trolley-load of a crane [1], [18] is given by where (x; z) are the coordinates of the load, m, which is connected to a trolley by a rope of length R and tension T . The trolley is at some distance D along the x-axis while the load is at an angle away from the vertical; see [1] and [18] for further details. As shown in these references, the system is differentially flat with a linearizing output given by y = [x; z] T .
The basic control problem is to carry the load m from (R1; D1) to (R 2 ; D 2 ) while minimizing oscillations at the end of the transport. Although the oscillations provide a natural way to formulate a cost functional, we use a slightly modified "indirect" approach suggested by The "high-level" control is obtained by [18] D( ) = y 1 ( ) 0 y1( )y2( ) y2( ) 0 g Recall that in our method, the derivatives are obtained by a simple matrix multiplication of the data at the LGL points. Fig. 1 displays plots in a form suitable for comparison with [1], where R1 = R2 = 10 m, D 1 = 0, D 2 = 20 m, and g = 9:8 m=s 2 , f = 10 s. An additional constraint z < g is also imposed to keep the rope tension positive. The number of LGL points were arbitrarily chosen to be N = 11. Although the shape of our plots is similar to that of Fliess et al., notice that our curves are a little different. This may be attributed to our use of higher-order polynomials where the extra degrees of freedom are used for optimization; hence, our method generates fewer oscillations as apparent from the plot of the vertical deviation angle, .

V. CONCLUSION AND FURTHER WORK
Pseudospectral (PS) methods offer a natural way to solve nonlinear control problems where the dynamics are described in terms of a differential-algebraic state space model. For flat systems, the optimal motion planning problem can be readily solved using PS methods. The computational ease derives from the use of higher order differentiation matrices and quadrature rules which transform the problem to a nonlinear programming problem with the values of the output variables at the quadrature nodes as the unknowns. While the notion of flatness is a promising idea, it is unclear at this stage whether optimal trajectories should be computed in (the flat) output space. In state-space, the boundary conditions are typically stated simply (e.g., linear boundary conditions) and have physical meaning. The flat output transforms these conditions to a possibly complex (e.g., nonlinear) set of end point conditions [compare (3) to (8)]. The same arguments hold for the transformation of the cost functional. Thus, it is possible that flatness parametrization might actually worsen real-time trajectory optimization. These issues are further elaborated in [19] with additional examples.