#  On the Rate of Convergence for the Pseudospectral Optimal Control of
Feedback Linearizable Systems ††thanks: The research was supported in part by
AFOSR and AFRL

Wei Kang
Department of Applied Mathematics
Naval Postgraduate School
Monterey, CA 93943
<EMAIL_ADDRESS>
###### Abstract

Over the last decade, pseudospectral (PS) computational methods for nonlinear
constrained optimal control have been applied to many industrial-strength
problems, notably the recent zero-propellant-maneuvering of the International
Space Station performed by NASA. In this paper, we prove a theorem on the rate
of convergence for the optimal cost computed using PS methods. It is a first
proved convergence rate in the literature of PS optimal control. In addition
to the high-order convergence rate, two theorems are proved for the existence
and convergence of the approximate solutions. This paper contains several
essential differences from existing papers on PS optimal control as well as
some other direct computational methods. The proofs do not use necessary
conditions of optimal control. Furthermore, we do not make coercivity type of
assumptions. As a result, the theory does not require the local uniqueness of
optimal solutions. In addition, a restrictive assumption on the cluster points
of discrete solutions made in existing convergence theorems are removed.

##  1 Introduction

Despite the fact that optimal control is one of the oldest problems in the
history of control theory, practical tools of solving nonlinear optimal
control problems are limited. Preferably, a feedback control law is derived
from a solution to the famously difficult Hamilton-Jacobi-Bellman (HJB)
equation. However, analytic solutions of this partial differential equation
can rarely be found for systems with nonlinear dynamics. Numerical
approximation of such solutions suffers from the well-known curse of
dimensionality and it is still an open problem for systems with moderately
high dimension. A practical alternative is to compute one optimal trajectory
at a time so that the difficulty of solving HJB equations is circumvented.
Then this open-loop optimal control can be combined with an inner-loop
tracking controller; or it can be utilized as a core instrument in a real-time
feedback control architecture such as a moving horizon feedback. A critical
challenge in this approach is to develop reliable and efficient computational
methods that generate the required optimal trajectories. In this paper, we
focus on some fundamental issues of pseudospectral computational optimal
control methods.

As a result of significant progress in large-scale computational algorithms
and nonlinear programming, the so-called direct computational methods have
become popular for solving nonlinear optimal control problems [1, 2, 17],
particularly in aerospace applications [16, 18]. In simple terms, in a direct
method, the continuous-time problem of optimal control is discretized, and the
resulting discretized optimization problem is solved by nonlinear programming
algorithms. Over the last decade, pseudospectral (PS) methods have emerged as
a popular direct methods for optimal control. They have been applied to many
industrial-strength problems, notably the recent attitude maneuvers of the
International Space Station performed by NASA. By following an attitude
trajectory developed using PS optimal control, the International Space Station
(ISS) was maneuvered 180 degrees on March 3, 2007, by using the gyroscopes
equipped on the ISS without propellant consumption. This single maneuver have
saved NASA about one million dollars’ worth of fuel [12]. The Legendre PS
optimal control method has already been developed into software named DIDO, a
MATLAB based package commercially available [19]. In addition, the next
generation of the OTIS software package [15] will have the Legendre PS method
as a problem solving option.

PS methods have been widely applied in scientific computation for models
governed by partial differential equations. The method is well known for being
very efficient in approximating solutions of differential equations. However,
despite its success and several decades of development, the intersection
between PS methods and nonlinear optimal control becomes an active research
area only after the mid-1990’s ([5, 6]). As yet, many fundamental theoretical
issues are still widely open. For the last decade, active research has been
carried out in the effort of developing a theoretical foundation for PS
optimal control methods. Among the research focuses, there are three
fundamental issues, namely the state and costate approximation, the existence
and convergence of approximate solutions, and the convergence rate. The
general importance of these issues is not limited to PS methods. They are
essential to other computational methods suchlike those based on Euler [10]
and Runge-Kutta [9] discretization. Similar to other direct computational
optimal control methods, PS method are based upon the Karush-Kuhn-Tucker (KKT)
conditions rather than the Pontryagin’s Minimum Principle (MPM). In [6] and
[8], a covector mapping was derived between the costate from KKT condition and
the costate from PMP. The covector mapping facilitates a verification and
validation of the computed solution. For the problem of convergence, some
theorems were published in [7]; and then the results were generalized in [13]
to problems with non-smooth control.

Among the three fundamental issues mentioned above, the most belated activity
of research is on the rate of convergence. In fact, there have been no results
published on the convergence rate for PS optimal control methods. Although
some results on the issue of convergence were proved in [7] and [13], a main
drawback of these results is the strong assumption in which the derivatives of
the discrete approximate solutions are required to converge uniformly. In this
paper, we prove a rate of convergence for the approximate optimal cost
computed using PS methods. Then, we prove theorems on existence and
convergence without the restrictive assumption made in [7] and [13]. In
addition to the high-order convergence rate addressed in Section 3, which is
the first proved convergence rate in the literature of PS optimal control,
this paper contains several essential differences from existing papers on PS
optimal control as well as some other direct computational methods. First of
all, the proof is not based on necessary conditions of optimal control.
Furthermore, we do not make coercivity type of assumptions. As a result, the
theory does not require the local uniqueness of optimal solutions. Therefore,
it is applicable to problems with multiple optimal solutions. Secondly, the
proof is not build on the bases of consistent approximation theory [17]. Thus,
we can remove the assumption in [7] and [13] on the existence of cluster
points for the derivatives of discrete solutions. The key that makes these
differences possible is that we introduce a set of sophisticated
regularization conditions in the discretization so that the computational
algorithm has a greater control of the boundedness of the approximate
solutions and their derivatives. Different from the existing results in the
literature of direct methods for optimal control, the desired boundedness is
achieved not by making assumptions on the original system, but by implementing
specially designed search region for the discrete problem of nonlinear
programming. This new boundary of search region automatically excludes
possible bad solutions that are numerically unstable.

The paper is organized as follows. In Section 2, the formulations of the
optimal control problem and its PS discretization are introduced. In 3, we
prove two theorems on the rate of convergence. In Section 4, two theorems on
the existence and convergence are proved.

##  2 Problem Formulation

For the rate of convergence, we focus on the following Bolza problem of
control systems in the feedback linearizable normal form. A more complicated
problem with constraints is studied in Section 4 for the existence and
convergence of approximate solutions.

Problem B: Determine the state-control function pair $(x(t),u(t))$,
$x\in\Re^{r}$ and $u\in\Re$, that minimizes the cost function

$\displaystyle J(x(\cdot),u(\cdot))$ $\displaystyle=$
$\displaystyle\int_{-1}^{1}F(x(t),u(t))\ dt+E(x(-1),x(1))$ (2.1)

subject to the following differential equations and initial condition

$\displaystyle\left\\{\begin{array}[]{lll}\dot{x}_{1}=x_{2}\\\ \;\;\;\vdots\\\
\dot{x}_{r-1}=x_{r}\\\ \dot{x}_{r}=f(x)+g(x)u\end{array}\right.$ (2.6)
$\displaystyle x(-1)=x_{0}$ (2.7)

where $x\in\Re^{r}$, $u\in\Re$, and $F:\Re^{r}\times\Re\to\Re$,
$E:\Re^{r}\times\Re^{r}\to\Re$, $f:\Re^{r}\to\Re$, and $g:\Re^{r}\to\Re$ are
all Lipschitz continuous functions with respect to their arguments. In
addition, we assume $g(x)\neq 0$ for all $x$.

Throughout the paper we make extensive use of Sobolev spaces, $W^{m,p}$, that
consists of functions, $\xi:[-1,1]\to\mathbb{R}$ whose $j$-th order weak
derivative, $\xi^{(j)}$, lies in $L^{p}$ for all $0\leq j\leq m$ with the
norm,

$\parallel\xi\parallel_{W^{m,p}}\quad=\sum_{j=0}^{m}\parallel\xi^{(j)}\parallel_{L^{p}}$

In this paper, we only consider the problems that have at least one optimal
solution in which $x_{r}^{\ast}(t)$ has bounded $m$-th order weak derivative,
i.e. $x_{r}^{\ast}(t)$ is in $W^{m,\infty}$. For some results, we assume
$m\geq 3$. For others, $m$ is smaller. Unless the term ‘strong derivative’ is
emphasized, all derivatives in the paper are in the weak sense.

The PS optimal control method addressed in this paper is an efficient direct
method. In typical direct methods, the original optimal control problem, not
the associated necessary conditions, is discretized to formulate a nonlinear
programming problem. The accuracy of the discretization is largely determined
by the accuracy of the underlying approximation method. Given any function
$f(t):[a,b]\rightarrow\Re$, a conventional method of approximation is to
interpolate at uniformly spaced nodes: $t_{0}=a$, $t_{1}=(b-a)/N$, $\cdots$,
$t_{N}=b$. However, it is known that uniform spacing is not efficient. More
sophisticated node selection methods are able to achieve significantly
improved accuracy with fewer nodes. It is important to emphasize that, for
optimal control problems, the rate of convergence is not merely an issue of
efficiency; more importantly it is about feasibility. An increased number of
nodes in discretization results in a higher dimension in the nonlinear
programming problem. A computational method becomes practically infeasible
when the dimension and complexity of the nonlinear programming exceed the
available computational power. In a PS approximation based on Legendre-Gauss-
Lobatto (LGL) quadrature nodes, a function $f(t)$ is approximated by $N$-th
order Lagrange polynomials using the interpolation at these nodes. The LGL
nodes, $t_{0}=-1<t_{1}<\cdots<t_{N}=1$, are defined by

$\begin{array}[]{llll}t_{0}=-1,\;\;t_{N}=1,\mbox{ and }\\\ \mbox{for
}k=1,2,\ldots,N-1,t_{k}\mbox{ are the roots of }\dot{L}_{N}(t)\end{array}$

where $\dot{L}_{N}(t)$ is the derivative of the $N$-th order Legendre
polynomial $L_{N}(t)$. The discretization works in the interval of $[-1,1]$.
An example of LGL nodes with $N=16$ is shown in Figure 1.

Figure 1: LGL nodes $N=16$

It was proved in approximation theory that the polynomial interpolation at the
LGL nodes converges to $f(t)$ under $L^{2}$ norm at the rate of $1/N^{m}$,
where $m$ is the smoothness of $f(t)$ (see for instance [4] Section 5.4). If
$f(t)$ is $C^{\infty}$, then the polynomial interpolation at the LGL nodes
converges at a spectral rate, i.e. it is faster than any given polynomial
rate. This is a very impressive convergence rate.

PS methods have been widely applied in scientific computation for models
governed by partial differential equations, such as complex fluid dynamics.
However, PS optimal control has several fundamental differences from the
computation of PDEs. Solving optimal control problems asks for the
approximation of several objects collectively, including the differential
equation that defines the control system, the integration in the cost
function, and the state and control trajectories. In addition to the various
types of approximations, a nonlinear programming must be applied to the
overall discretized optimization problem to find an approximate optimal
control. All these factors may deteriorate the final approximate solution. The
existing theory of PS approximation of differential equations is not
applicable. New theory needs to be developed for the existence, convergence,
and the rate of convergence for optimal control problems.

In the following, we introduce the notations used in this paper. Then, the
discretized nonlinear programming problem is formulated. In a PS optimal
control method, the state and control functions, $x(t)$ and $u(t)$, are
approximated by $N$-th order Lagrange polynomials based on the interpolation
at the LGL quadrature nodes. In the discretization, the state variables are
approximated by the vectors $\bar{x}^{Nk}\in\Re^{r}$, i.e.

$\bar{x}^{Nk}=\left[\begin{array}[]{cccccccccccccc}\bar{x}_{1}^{Nk}\\\
\bar{x}_{2}^{Nk}\\\ \vdots\\\ \bar{x}_{r}^{Nk}\end{array}\right]$

is an approximation of $x(t_{k})$. Similarly, $\bar{u}^{Nk}$ is the
approximation of $u(t_{k})$. Thus, a discrete approximation of the function
$x_{i}(t)$ is the vector

$\bar{x}_{i}^{N}=\left[\begin{array}[]{cccccccccccccc}\bar{x}_{i}^{N1}&\bar{x}_{i}^{N2}&\cdots&\bar{x}_{i}^{NN}\end{array}\right]$

A continuous approximation is defined by its polynomial interpolation, denoted
by $x_{i}^{N}(t)$, i.e.

$\displaystyle x_{i}(t)$ $\displaystyle\approx$ $\displaystyle
x_{i}^{N}(t)=\sum_{k=0}^{N}\bar{x}_{i}^{Nk}\phi_{k}(t),$ (2.8)

where $\phi_{k}(t)$ is the Lagrange interpolating polynomial [4]. Instead of
polynomial interpolation, the control input is approximated by the following
non-polynomial interpolation

$\displaystyle
u^{N}(t)=\displaystyle\frac{\dot{x}_{r}^{N}(t)-f(x^{N}(t))}{g(x^{N}(t))}$
(2.9)

In the notations, the discrete variables are denoted by letters with an upper
bar, such as $\bar{x}^{Nk}_{i}$ and $\bar{u}^{Nk}$. If $k$ in the superscript
and/or $i$ in the subscript are missing, it represents the corresponding
vector or matrix in which the indices run from minimum to maximum. For
example,

$\displaystyle\bar{x}^{N}_{i}$ $\displaystyle=$
$\displaystyle\left[\begin{array}[]{cccccccccccccc}\bar{x}_{i}^{N0}&\bar{x}_{i}^{N1}&\cdots&\bar{x}_{i}^{NN}\end{array}\right]$
$\displaystyle\bar{x}^{Nk}$ $\displaystyle=$
$\displaystyle\left[\begin{array}[]{cccccccccccccc}\bar{x}_{1}^{Nk}\\\
\bar{x}_{2}^{Nk}\\\ \vdots\\\ \bar{x}_{r}^{Nk}\end{array}\right]$
$\displaystyle\bar{x}^{N}$ $\displaystyle=$
$\displaystyle\left[\begin{array}[]{cccccccccccccc}\bar{x}_{1}^{N0}&\bar{x}_{1}^{N1}&\cdots&\bar{x}_{1}^{NN}\\\
\bar{x}_{2}^{N0}&\bar{x}_{2}^{N1}&\cdots&\bar{x}_{2}^{NN}\\\
\vdots&\vdots&\vdots&\vdots\\\
\bar{x}_{r}^{N0}&\bar{x}_{r}^{N1}&\cdots&\bar{x}_{r}^{NN}\end{array}\right]$

Similarly,

$\bar{u}^{N}=\left[\begin{array}[]{cccccccccccccc}\bar{u}^{N0}&\bar{u}^{N1}&\cdots&\bar{u}^{NN}\end{array}\right]$

Given a discrete approximation of a continuous function, the interpolation is
denoted by the same notation without the upper bar. For example,
$x_{i}^{N}(t)$ in (2.8), $u^{N}(t)$ in (2.9). The superscript $N$ represents
the number of LGL nodes used in the approximation. Throughout this paper, the
interpolation of $(\bar{x}^{N},\bar{u}^{N})$ is defined by (2.8)-(2.9), in
which $u^{N}(t)$ is not necessarily a polynomial. It is proved in Lemma 5 that
(2.9) is indeed an interpolation.

Existing results in the analysis of spectral methods show that PS method is an
approach that is easy and accurate in the approximation of smooth functions,
integrations, and differentiations, all critical to optimal control problems.
For differentiation, the derivative of $x^{N}_{i}(t)$ at the LGL node $t_{k}$
is easily computed by the following matrix multiplication [4]

$\displaystyle\left[\begin{array}[]{cccccccccccccc}\dot{x}_{i}^{N}(t_{0})&\dot{x}_{i}^{N}(t_{1})&\cdots&\dot{x}_{i}^{N}(t_{N})\end{array}\right]^{T}=D(\bar{x}^{N}_{i})^{T}$
(2.14)

where the $(N+1)\times(N+1)$ differentiation matrix $D$ is defined by

$\displaystyle D_{ik}$ $\displaystyle=$
$\displaystyle\left\\{\begin{array}[]{ll}\frac{L_{N}(t_{i})}{L_{N}(t_{k})}\frac{1}{t_{i}-t_{k}},&\mbox{if}\
\ i\neq k;\\\ \\\ -\frac{N(N+1)}{4},&\mbox{if}\ \ i=k=0;\\\ \\\
\frac{N(N+1)}{4},&\mbox{if}\ \ i=k=N;\\\ \\\
0,&\mbox{otherwise}\end{array}\right.$

The cost functional $J[x(\cdot),u(\cdot)]$ is approximated by the Gauss-
Lobatto integration rule,

$\displaystyle J[x(\cdot),u(\cdot)]\ \approx\
\bar{J}^{N}(\bar{x}^{N},\bar{u}^{N})$ $\displaystyle=$
$\displaystyle\sum_{k=0}^{N}F(\bar{x}^{Nk},\bar{u}^{Nk})w_{k}+E(\bar{x}^{N0},\bar{x}^{NN})$

where $w_{k}$ are the LGL weights defined by

$\displaystyle w_{k}$ $\displaystyle=$
$\displaystyle\frac{2}{N(N+1)}\frac{1}{[L_{N}(t_{k})]^{2}},$

The approximation is so accurate that it has zero error if the integrand
function is a polynomial of degree less than or equal to $2N-1$, a degree that
is almost a double of the number of nodes [4]. Now, we are ready to define
Problem ${\rm B}^{\rm N}$, a PS discretization of Problem B.

For any integer $m_{1}>0$, let
$\\{a_{0}^{N}(m_{1}),a_{1}^{N}(m_{1}),\cdots,a_{N-r-m_{1}+1}^{N}(m_{1})\\}$
denote the coefficients in the Legendre polynomial expansion for the
interpolation polynomial of the vector $\bar{x}_{r}^{N}(D^{T})^{m_{1}}$. Note
that the interpolation of $\bar{x}_{r}^{N}(D^{T})^{m_{1}}$ equals the
polynomial of $\frac{d^{m_{1}}x_{r}^{N}(t)}{dt^{m_{1}}}$. Thus, there are only
$N-r-m_{1}+2$ nonzero spectral coefficients because it is proved in Section 3
that the order of $\frac{d^{m_{1}}x^{N}_{r}(t)}{dt^{m_{1}}}$ is at most degree
of $N-r-m_{1}+1$. These coefficients depend linearly on $\bar{x}_{r}^{N}$ [3],

$\begin{array}[]{llllllllll}{\tiny\left[\begin{array}[]{cccccccccccccc}a^{N}_{0}(m_{1})\\\
\vdots\\\ a^{N}_{N-r-m_{1}+1}(m_{1})\end{array}\right]=}\\\
{\tiny\left[\begin{array}[]{cccccccccccccc}\frac{1}{2}&&\\\ &\ddots&\\\ &&N-r-
m_{1}+1+\frac{1}{2}\end{array}\right]\left[\begin{array}[]{cccccccccccccc}L_{0}(t_{0})&\cdots&L_{0}(t_{N})\\\
&\vdots&\\\ L_{N-r-m_{1}+1}(t_{0})&\cdots&L_{N-r-
m_{1}+1}(t_{N})\end{array}\right]\left[\begin{array}[]{cccccccccccccc}w_{0}&&\\\
&\ddots&\\\ &&w_{N}\\\
\end{array}\right]D^{m_{1}}\left[\begin{array}[]{cccccccccccccc}\bar{x}^{N0}_{r}\\\
\vdots\\\ \bar{x}^{NN}_{r}\\\ \end{array}\right]}\end{array}$ (2.16)

The PS discretization of Problem ${\rm B}^{\rm N}$is defined as follows.

Problem ${\bf B}^{\bf N}$: Find $\bar{x}^{Nk}\in\Re^{r}$ and
$\bar{u}^{Nk}\in\Re$, $k\ =\ 0,1,\ldots,N$, that minimize

$\displaystyle\bar{J}^{N}(\bar{x}^{N},\bar{u}^{N})$ $\displaystyle=$
$\displaystyle\sum_{k=0}^{N}F(\bar{x}^{Nk},\bar{u}^{Nk})w_{k}+E(\bar{x}^{N0},\bar{x}^{NN})$
(2.17)

subject to

$\displaystyle\left\\{\begin{array}[]{rcl}D(\bar{x}_{1}^{N})^{T}&=&(\bar{x}_{2}^{N})^{T}\\\
D(\bar{x}_{2}^{N})^{T}&=&(\bar{x}_{3}^{N})^{T}\\\ &\vdots&\\\
D(\bar{x}_{r-1}^{N})^{T}&=&(\bar{x}_{r}^{N})^{T}\\\
D(\bar{x}_{r}^{N})^{T}&=&\left[\begin{array}[]{cccccccccccccc}f(\bar{x}^{N0})+g(\bar{x}^{N0})\bar{u}^{N0}\\\
\vdots\\\ f(\bar{x}^{NN})+g(\bar{x}^{NN})\bar{u}^{NN}\end{array}\right]\\\
\end{array}\right.$ (2.26) $\displaystyle\bar{x}^{N0}=x_{0}$ (2.27)
$\displaystyle\underline{{\boldsymbol{b}}}\leq\left[\begin{array}[]{cccccccccccccc}\bar{x}^{Nk}\\\
\bar{u}^{Nk}\end{array}\right]\ \leq\ \bar{\boldsymbol{b}},\;\;\;\;\mbox{ for
all }0\leq k\leq N$ (2.30)
$\displaystyle\underline{{\boldsymbol{b}}}_{j}\leq\left[\begin{array}[]{cccccccccccccc}1&0&\cdots&0\end{array}\right]D^{j}(\bar{x}_{r}^{N})^{T}\
\leq\bar{\boldsymbol{b}}_{j},\mbox{ if }1\leq j\leq m_{1}-1\mbox{ and
}m_{1}\geq 2$ (2.32) $\displaystyle\displaystyle\sum_{n=0}^{N-r-
m_{1}+1}|a^{N}_{n}(m_{1})|\leq{\boldsymbol{d}}$ (2.33)

Comparing to Problem B, (2.26) is the discretization of the control system
defined by the differential equation. The regularization condition (2.33)
assures that the derivative of the interpolation up to the order of $m_{1}$ is
bounded. It is proved in the following sections that the integer $m_{1}$ is
closely related to the convergence rate. The inequalities (2.30), (2.32) and
(2.33) are regularization conditions that do not exist in Problem B. It is
proved in the next few sections that these additional constraints do not
affect the feasibility of Problem ${\rm B}^{\rm N}$. Therefore, it does not
put an extra limit to the family of problems to be solved.

In searching for a discrete optimal solution, it is standard for software
packages of nonlinear programming to require a search region. Typically, the
search region is defined by a constraint (2.30). However, this box-shaped
region may contain solutions that are not good approximations of the
continuous-time solution. To guarantee the rate of convergence, the search
region is refined to a smaller one by imposing constraints (2.32) and (2.33).
It is proved in this paper that there always exist feasible solutions that
satisfy all the constraints and the optimal cost converges, provided the upper
and lower bounds are large enough. In (2.32),
$\underline{{\boldsymbol{b}}}_{j}$ and $\bar{\boldsymbol{b}}_{j}$ represent
the bounds of initial derivatives. In (2.33), ${\boldsymbol{d}}$ is the bound
determined by $x_{r}^{(m_{1}+1)}$ satisfying the inequality (3.6). Without
known the optimal solution, these bounds of search region have to be estimated
before computation or they are determined by numerical experimentations. The
constraints (2.32) and (2.33) are necessary to avoid the restrictive
consistent approximation assumption made in [7]. At a more fundamental level,
the order of derivatives, $m_{1}$ in (2.33), determines the convergence rate
of the approximate optimal control. Another interesting fact that amply
justify these additional constraints is that Problem ${\rm B}^{\rm N}$ may not
even have an optimal solution if we do not enforce (2.30). This is shown by
the following counter example.

###### Example 1

Consider the following problem of optimal control.

$\displaystyle\min_{(x(\cdot),u(\cdot))}\int_{-1}^{1}\displaystyle\frac{(x(t)-u(t))^{2}}{u(t)^{4}}dt$
$\displaystyle\dot{x}=u$ (2.34) $\displaystyle x(-1)=e^{-1}$

It is easy to check that the optimal solution is

$\displaystyle u=e^{t},$ $\displaystyle x(t)=e^{t}$

and the optimal cost value is zero. Although the solution to the problem
(2.34) is simple and analytic, the PS discretization of (2.34) does not have
an optimal solution if the constraint (2.30) is not enforced. To prove this
claim, consider the PS discretization,

$\displaystyle\min_{(\bar{x}^{N},\bar{u}^{N})}\bar{J}^{N}(\bar{x}^{N},\bar{u}^{N})=\displaystyle\sum_{k=0}^{N}\displaystyle\frac{(\bar{x}^{Nk}-D_{k}(\bar{x}^{N})^{T})^{2}}{\left(D_{k}(\bar{x}^{N})^{T}\right)^{4}}w_{k}$
$\displaystyle D(\bar{x}^{N})^{T}=(\bar{u}^{N})^{T}$ (2.35)
$\displaystyle\bar{x}^{N0}=e^{-1}$

where $D_{k}$ is the $k$th row of the differentiation matrix $D$. Let
$x^{N}(t)$ be the interpolation polynomial of $\bar{x}^{N}$, then it is
obvious that

$x^{N}(t)-\dot{x}^{N}(t)\not\equiv 0$

Thus, there exists $k$ so that

$\bar{x}^{Nk}-D_{k}(\bar{x}^{N})^{T}\neq 0$

So,

$\begin{array}[]{llllllllll}\bar{J}^{N}(\bar{x}^{N},\bar{u}^{N})>0\end{array}$
(2.36)

for all feasible pairs $(\bar{x}^{N},\bar{u}^{N})$. For any $\alpha>0$, define

$\bar{x}^{Nk}=e^{-1}+\alpha(t_{k}+1)$

The interpolation of $\bar{x}^{N}$ is the linear polynomial

$x^{N}(t)=e^{-1}+\alpha(t+1)$

Then,

$D_{k}(\bar{x}^{N})^{T}=\dot{x}^{N}(t_{k})=\alpha$

The cost function is

$\displaystyle\bar{J}^{N}(\bar{x}^{N},\bar{u}^{N})$ $\displaystyle=$
$\displaystyle\displaystyle\sum_{k=0}^{N}\displaystyle\frac{(e^{-1}+\alpha(t_{k}+1)-\alpha)^{2}}{\alpha^{4}}w_{k}$
$\displaystyle=$
$\displaystyle\displaystyle\sum_{k=0}^{N}\displaystyle\frac{(e^{-1}+\alpha
t_{k})^{2}}{\alpha^{4}}w_{k}$ $\displaystyle\leq$
$\displaystyle\displaystyle\sum_{k=0}^{N}\displaystyle\frac{(e^{-1}+\alpha)^{2}}{\alpha^{4}}w_{k}$
$\displaystyle=$ $\displaystyle
2\displaystyle\frac{(e^{-1}+\alpha)^{2}}{\alpha^{4}}$

Therefore, $\bar{J}^{N}(\bar{x}^{N},\bar{u}^{N})$ can be arbitrarily small as
$\alpha$ approaches $\infty$. However, $\bar{J}^{N}(\bar{x}^{N},\bar{u}^{N})$
is always positive as shown by (2.36). We conclude that the discretization
(2.35) has no minimum value for $\bar{J}^{N}(\bar{x}^{N},\bar{u}^{N})$.

##  3 Convergence Rate

Given a solution to Problem ${\rm B}^{\rm N}$, we use (2.9) to approximate the
optimal control. In this section we prove that, under this approximate optimal
control, the value of the cost function converges to the optimal cost of
Problem B as the number of nodes is increased. More importantly, we can prove
a high-order rate of convergence. In the literature, it has been proved that
PS methods have a spectral rate when approximating $C^{\infty}$ functions,
i.e. the rate is faster than any polynomial rate. However, there are no
results in the literature thus far on the convergence rate of PS optimal
control. Meanwhile, in many problems solved by PS optimal control we clearly
observed a rate of high-order in the convergence. In this section, we prove a
convergence rate that depends on the smoothness of the optimal control. More
specifically, the rate is about $\frac{1}{N^{2m/3-1}}$, where $m$ is defined
by the smoothness of the optimal trajectory. If the cost function can be
accurately computed, then the convergence rate is improved to
$\frac{1}{N^{2m-1}}$. In the special case of $C^{\infty}$, it is proved that
PS method is able to converge faster than any given polynomial rate. Before we
introduce the main theorems of this section, the following example in [7] is
briefly presented to show the rapid convergence of the PS optimal control
method.

###### Example 2

Consider the following nonlinear optimal control problem:

$\displaystyle\left\\{\begin{array}[]{lrl}{\rm
Minimize}&J[x(\cdot),u(\cdot)]=&4x_{1}(2)+x_{2}(2)+4\displaystyle\int_{0}^{2}u^{2}(t)\
dt\\\ {\rm Subject\ to}&\dot{x}_{1}(t)=&x_{2}^{3}(t)\\\
&\dot{x}_{2}(t)=&u(t)\\\ &(x_{1}(0),x_{2}(0))=&(0,1)\end{array}\right.$ (3.5)

It can be shown that the exact optimal control is defined by
$u^{\ast}(t)=-\frac{8}{(2+t)^{3}}$. For this problem, the PS method achieves
the accuracy in the magnitude of $10^{-8}$ with only 18 nodes. A detail
comparison of the PS method with some other discretization methods are
addressed in [7]. From Figure 2 in logarithmically scaled coordinates, it is
obvious that the computation using the PS method converges exponentially.

Figure 2: Error vs number of the nodes for the pseudospectral method

Problem ${\rm B}^{\rm N}$ has several bounds in its definition,
$\underline{{\boldsymbol{b}}}$, $\bar{\boldsymbol{b}}$,
$\underline{{\boldsymbol{b}}}_{j}$, $\bar{\boldsymbol{b}}_{j}$, and
${\boldsymbol{d}}$. These bounds can be selected from a range determined by
Problem B. The constraints $\underline{{\boldsymbol{b}}}$ and
$\bar{\boldsymbol{b}}$ are lower and upper bounds so that the optimal
trajectory of Problem B is contained in the interior of the region. Suppose
Problem B has an optimal solution $(x^{\ast}(t),u^{\ast}(t))$ in which
$(x_{r}^{\ast}(t))^{(m)}$ has bounded variation for some $m\geq 3$, where
$x_{r}^{\ast}(t)$ is the $r$th component of the optimal trajectory. Suppose
$m_{1}$ in Problem ${\rm B}^{\rm N}$ satisfies $2\leq m_{1}\leq m-1$. Then, we
can select the bounds $\underline{{\boldsymbol{b}}}_{j}$ and
$\bar{\boldsymbol{b}}_{j}$ so that $(x_{r}^{\ast}(t))^{(j)}$ is contained in
the interior of the region. For ${\boldsymbol{d}}$, we assume

$\displaystyle{\boldsymbol{d}}$ $\displaystyle>$
$\displaystyle\displaystyle\frac{6}{\sqrt{\pi}}(U(\left.x_{r}^{\ast}\right.^{(m_{1}+1)})+V(\left.x_{r}^{\ast}\right.^{(m_{1}+1)}))\zeta(3/2)$
(3.6)

where $U(\left.x_{r}^{\ast}\right.^{(m_{1}+1)})$ is the upper bound and
$V(\left.x_{r}^{\ast}\right.^{(m_{1}+1)})$ is the total variation of
$\left.x_{r}^{\ast}\right.^{(m_{1}+1)}(t)$; and $\zeta(s)$ is the $\zeta$
function defined by

$\displaystyle\zeta(s)=\sum_{k=1}^{\infty}\displaystyle\frac{1}{k^{s}}$ (3.7)

If all the bounds are selected as above, then it is proved in Section 4 that
Problem ${\rm B}^{\rm N}$ is always feasible provided $m\geq 2$. Note that in
practical computation, $\underline{{\boldsymbol{b}}}$, $\bar{\boldsymbol{b}}$,
$\underline{{\boldsymbol{b}}}_{j}$, $\bar{\boldsymbol{b}}_{j}$, and
${\boldsymbol{d}}$ are unknown. They must be estimated based upon experience
or other information about the system.

###### Theorem 1

Suppose Problem B has an optimal solution $(x^{\ast}(t),u^{\ast}(t))$ in which
the $m$-th order derivative $(x_{r}^{\ast}(t))^{(m)}$ has a bounded variation
for some $m\geq 3$. In Problem ${\rm B}^{\rm N}$, select $m_{1}$ and $\alpha$
so that $1\leq m_{1}\leq m-1$ and $0<\alpha<m_{1}-1$. Suppose $f(\cdot)$,
$g(\cdot)$, $F(\cdot)$, and $E(x_{0},\cdot)$ are $C^{m}$ and globally
Lipschitz. Suppose all other bounds in Problem ${\rm B}^{\rm N}$ are large
enough. Given any sequence

$\begin{array}[]{llllllllll}\\{(\bar{x}^{\ast N},\bar{u}^{\ast N})\\}_{N\geq
N_{1}}\end{array}$ (3.8)

of optimal solutions of Problem ${\rm B}^{\rm N}$. Then the approximate cost
converge to the optimal value at the following rate

$\displaystyle\left|J(x^{\ast}(\cdot),u^{\ast}(\cdot))-J(x^{\ast
N}(\cdot),u^{\ast N}(\cdot))\right|$ $\displaystyle\leq$
$\displaystyle\displaystyle\frac{M_{1}}{(N-r-
m_{1}-1)^{2m-2m_{1}-1}}+\displaystyle\frac{M_{2}}{N^{\alpha}}$ (3.9)
$\displaystyle\left|J(x^{\ast}(\cdot),u^{\ast}(\cdot))-\bar{J}^{N}(\bar{x}^{\ast
N},\bar{u}^{\ast N})\right|$ $\displaystyle\leq$
$\displaystyle\displaystyle\frac{M_{1}}{(N-r-
m_{1}-1)^{2m-2m_{1}-1}}+\displaystyle\frac{M_{2}}{N^{\alpha}}$ (3.10)

where $M_{1}$ and $M_{2}$ are some constants independent of $N$. In (3.9),
$(x^{\ast N}(t),u^{\ast N}(t))$ is the interpolation of (3.8) defined by
(2.8)-(2.9). In fact, $x^{\ast N}(t)$ is the trajectory of (2.6) under the
control input $u^{\ast N}(t)$.

Theorem 1 implies that the costs of any sequence of discrete optimal solutions
must converge to the optimal cost of Problem B, no matter the sequence of the
discrete state and control trajectories converge or not. In other words, it is
possible that the sequence of discrete optimal controls does not converge to a
unique continuous-time control; meanwhile the costs using these approximate
optimal controls converge to the true optimal cost of Problem B. Therefore,
this theorem does not require the local uniqueness of solutions for Problem B.
This is different from many existing convergence theorems of computational
optimal control, in which a unique optimal solution and coercivity are
assumed. This is made possible because the proofs in this paper do not rely on
the necessary conditions of optimal control. The key idea in the proof is to
shape the search region in Problem ${\rm B}^{\rm N}$by regulating the discrete
solutions using (2.30)-(2.32)-(2.33). We would like to emphasize that the
regulation constraints are added to the discretized problem, not the original
Problem B. So, the constraints do not restrict the problem to be solved, and
they do not need to be verified before computation. In addition, increasing
the number of constraints results in smaller search region for an optimal
solution.

###### Remark 3.1

If $f(\cdot)$, $g(\cdot)$, $F(\cdot)$ and $x_{r}^{\ast}(t)$ are $C^{\infty}$,
then we can select $m$ and $m_{1}$ arbitrarily large. In this case, we can
make the optimal cost of Problem ${\rm B}^{\rm N}$converge faster than any
given polynomial rate.

###### Remark 3.2

From (3.9) and (3.10), the convergence rate is determined by $m$, the
smoothness of the optimal trajectory, and $m_{1}$, the order in the regulation
of discrete solutions. While $m$ is a property of Problem B that cannot be
changed, $m_{1}$ in Problem ${\rm B}^{\rm N}$ can be selected within a range.
However, the errors in (3.9) and (3.10) have two parts, one is an increasing
function of $m_{1}$ and the other is a decreasing function of $m_{1}$. In
Corollary 1, we show an optimal selection of $m_{1}$ to maximize the combined
convergence rate.

The proof is convoluted involving results from several different areas,
including nonlinear functional analysis, orthogonal polynomials, and
approximation theory. First, we introduce the concept of Fréchet derivative.
Let us consider the continuous cost function, $J(x(\cdot),u(\cdot))$, subject
to (2.6)-(2.7) as a nonlinear functional of $u(\cdot)$, denoted by ${\cal
J}(u)$. For any $u$ in the Banach space $W^{m-1,\infty}$, suppose there exists
a linear bounded operator ${\cal L}$: $W^{m-1,\infty}\rightarrow\Re$ such that

$|{\cal J}(u+\Delta u)-{\cal J}(u)-L\Delta u|=o(||\Delta
u||_{W^{m-1,\infty}})$

for all $u+\Delta u$ in an open neighborhood of $u$ in $W^{m-1,\infty}$. Then,
${\cal L}$ is called the Fréchet derivative of ${\cal J}(u)$ at $u$, denoted
by ${\cal J}^{\prime}(u)=\cal L$. If ${\cal J}^{\prime}(u)$ exists at all
points in an open subset of $W^{m-1,\infty}$, then ${\cal J}^{\prime}(u)$ is a
functional from this open set to the Banach space $L(W^{m-1,\infty},\Re)$ of
all bounded linear operators. If this new functional has a Fréchet derivative,
then it is called the second order Fréchet derivative, denoted by ${\cal
J}^{\prime\prime}(u)$. The following lemma is standard in nonlinear functional
analysis [21].

###### Lemma 1

Suppose ${\cal J}$ takes a local minimum value at $u^{\ast}$. Suppose ${\cal
J}$ has second order Fréchet derivative at $u^{\ast}$. Then,

${\cal J}(u^{\ast}+\Delta u)=({\cal J}^{\prime\prime}(u^{\ast})\Delta u)\Delta
u+o(||\Delta u||^{2})$

The rate of convergence for the spectral coefficients can be estimated by the
following Jackson’s Theorem.

###### Lemma 2

(Jackson’s Theorem [20]) Let $h(t)$ be of bounded variation in $[-1,1]$.
Define

$H(t)=H(-1)+\displaystyle{\int_{-1}^{t}}h(s)ds$

then $\\{a_{n}\\}_{n=0}^{\infty}$, the sequence of spectral coefficients of
$H(t)$, satisfies the following inequality

$a_{n}<\displaystyle\frac{6}{\sqrt{\pi}}(U(h(t))+V(h(t)))\displaystyle\frac{1}{n^{3/2}}$

for $n\geq 1$.

Given a continuous function $h(t)$ defined on $[-1,1]$. Let $\hat{p}^{N}(t)$
be the best polynomial of degree $N$, i.e. the $N$th order polynomial with the
smallest distance to $h(t)$ under $||\cdot||_{\infty}$ norm. Let $I_{N}h(t)$
be the polynomial interpolation using the value of $h(t)$ at the LGL nodes.
Then, we have the following inequality from the theory of approximation and
orthogonal polynomials [4], [11].

###### Lemma 3

$||h(t)-I_{N}h||_{\infty}\leq(1+\Lambda_{N})||h(t)-\hat{p}^{N}(t)||_{\infty}$

where $\Lambda_{N}$ is called Lebesgue constant. It satisfies

$\Lambda_{N}\leq\displaystyle\frac{2}{\pi}log(N+1)+0.685\cdots$

The best polynomial approximation represents the closest polynomial to a
function under $||\cdot||_{\infty}$. The error can be estimated by the
following Lemma [4].

###### Lemma 4

(1) Suppose $h(t)\in W^{m,\infty}$. Let $\hat{p}^{N}(t)$ be the best
polynomial approximation. Then

$||\hat{p}^{N}(t)-h(t)||_{\infty}\leq\displaystyle\frac{C}{N^{m}}||h(t)||_{W^{m,\infty}}$

for some constant $C$ independent of $h(t)$, $m$ and $N$.

(2) If $h(t)\in W^{m,2}$, then

$||h(t)-P_{N}h(t)||_{\infty}\leq\displaystyle\frac{C||h(t)||_{W^{m,2}}}{N^{m-3/4}}$

where $P_{N}h$ is the N-th order truncation of the Legendre series of $h(t)$.
(3) If $h(t)$ has the $m$-th order strong derivative with a bounded variation,
then

$||h(t)-P_{N}h(t)||_{\infty}\leq\displaystyle\frac{CV(h^{(m)}(t))}{N^{m-1/2}}$

The following lemmas are proved specifically for PS optimal control methods.
Similar results can be found in [14] except that some assumptions on $m_{1}$
are relaxed.

###### Lemma 5

([14]) (i) For any trajectory, $(\bar{x}^{N},\bar{u}^{N})$, of the dynamics
(2.26), the pair $(x^{N}(t),u^{N}(t))$ defined by (2.8)-(2.9) satisfies the
differential equations defined in (2.6). Furthermore,

$\displaystyle\bar{x}^{Nk}=x^{N}(t_{k}),\;\bar{u}^{Nk}=u^{N}(t_{k}),\;\mbox{for
}k=0,1,\cdots,N$ (3.11)

(ii) For any pair $(x^{N}(t),u^{N}(t))$ in which $x^{N}(t)$ consists of
polynomials of degree less than or equal to $N$ and $u^{N}(t)$ is a function,
if $(x^{N}(t),u^{N}(t))$ satisfies the differential equations in (2.6), then
$(\bar{x}^{N},\bar{u}^{N})$ defined by (3.11) satisfies (2.26).

(iii) If $(\bar{x}^{N},\bar{u}^{N})$ satisfies (2.26), then the degree of
$x_{i}^{N}(t)$ is less than or equal to $N-i+1$.

Proof.  (i) Suppose $(\bar{x}^{N},\bar{u}^{N})$ satisfies the equations in
(2.26). Because $x^{N}(t)$ is the polynomial interpolation of $\bar{x}^{N}$,
and because of equations (2.14), we have

$\begin{array}[]{llllllllll}\left[\begin{array}[]{cccccccccccccc}\dot{x}_{i}^{N}(t_{0})&\dot{x}_{i}^{N}(t_{1})&\cdots&\dot{x}_{i}^{N}(t_{N})\end{array}\right]\\\
=\bar{x}_{i}^{N}D^{T}\\\ =\bar{x}_{i+1}^{N}\\\
=\left[\begin{array}[]{cccccccccccccc}x_{i+1}^{N}(t_{0})&x_{i+1}^{N}(t_{1})&\cdots&x_{i+1}^{N}(t_{N})\end{array}\right]\end{array}$

Therefore, the polynomials $\dot{x}^{N}_{i}(t)$ and $x^{N}_{i+1}(t)$ must
equal each other because they coincide at $N+1$ points and because the degrees
of $x_{i}^{N}(t)$ and $x_{i+1}^{N}(t)$ are less than or equal to $N$. In
addition, (2.9), the definition of $u^{N}(t)$, implies the last equation in
(2.6). So, the pair $(x^{N}(t),u^{N}(t))$ satisfies all equations in (2.6).
Now, we prove (3.11). Because $x^{N}(t)$ is an interpolation of $\bar{x}^{N}$,
we know $\bar{x}^{Nk}=x^{N}(t_{k})$ for $0\leq k\leq N$. From (2.9),

$\displaystyle u^{N}(t_{k})$ $\displaystyle=$
$\displaystyle\displaystyle\frac{\dot{x}_{r}^{N}(t_{k})-f(x^{N}(t_{k}))}{g(x^{N}(t_{k}))}$
(3.12) $\displaystyle=$
$\displaystyle\displaystyle\frac{\dot{x}^{N}_{r}(t_{k})-f(\bar{x}^{Nk})}{g(\bar{x}^{Nk})}$

Because of (2.14), we have

$\left[\begin{array}[]{cccccccccccccc}\dot{x}_{r}^{N}(t_{0})&\dot{x}_{r}^{N}(t_{1})&\cdots&\dot{x}_{r}^{N}(t_{N})\end{array}\right]^{T}=D(\bar{x}_{r}^{N})^{T}$

Therefore, (3.12) is equivalent to

$\displaystyle\left[\begin{array}[]{cccccccccccccc}u^{N}(t_{0})&u^{N}(t_{1})&\cdots&u^{N}(t_{N})\end{array}\right]^{T}$
$\displaystyle=$
$\displaystyle\mbox{diag}\left(\displaystyle\frac{1}{g(\bar{x}^{N0})},\cdots,\displaystyle\frac{1}{g(\bar{x}^{NN})}\right)\left(D(\bar{x}_{r}^{N})^{T}-\left[\begin{array}[]{cccccccccccccc}f(\bar{x}^{N0})\\\
\vdots\\\ f(\bar{x}^{NN})\end{array}\right]\right)$

Comparing to the last equation in (2.26), it is obvious that
$u^{N}(t_{k})=\bar{u}^{Nk}$. So, (3.11) holds true. Part (i) is proved.

(ii) Assume $(x^{N}(t),u^{N}(t))$ satisfies the differential equations in
(2.6). Because $x^{N}(t)$ are polynomials, (2.14) implies

$\displaystyle\bar{x}_{i}^{N}D^{T}$
$\displaystyle=\left[\begin{array}[]{cccccccccccccc}\dot{x}_{i}^{N}(t_{0})&\dot{x}_{i}^{N}(t_{1})&\cdots&\dot{x}_{i}^{N}(t_{N})\end{array}\right]$
(3.15)
$\displaystyle=\left[\begin{array}[]{cccccccccccccc}x_{i+1}^{N}(t_{0})&x_{i+1}^{N}(t_{1})&\cdots&x_{i+1}^{N}(t_{N})\end{array}\right]$
(3.17) $\displaystyle=\bar{x}_{i+1}^{N}$

Furthermore,

$\displaystyle\bar{x}^{N}_{r}D^{T}$ $\displaystyle=$
$\displaystyle\left[\begin{array}[]{cccccccccccccc}\dot{x}_{r}^{N}(t_{0})&\dot{x}_{r}^{N}(t_{1})&\cdots&\dot{x}_{r}^{N}(t_{N})\end{array}\right]$
(3.19) $\displaystyle=$
$\displaystyle\left[\begin{array}[]{cccccccccccccc}f(x^{N}(t_{0}))+g(x^{N}(t_{0}))u^{N}(t_{0})&\cdots&f(x^{N}(t_{N}))+g(x^{N}(t_{N}))u^{N}(t_{N})\end{array}\right]$
(3.21)

Equations (3.15) and (3.19) imply that $(\bar{x}^{N},\bar{u}^{N})$ satisfies
(2.26). Part (ii) is proved.

(iii) We know that the degree of $x_{1}^{N}(t)$, the interpolation polynomial,
is less than or equal to $N$. From (i), we know
$x_{2}^{N}(t)=\dot{x}_{1}^{N}(t)$. Therefore, the degree of $x_{2}^{N}(t)$
must be less than or equal to $N-1$. In general, the degree of $x_{i}^{N}(t)$
is less than or equal to $N-i+1$. $\Box$

###### Lemma 6

Suppose $\\{(\bar{x}^{N},\bar{u}^{N})\\}_{N=N_{1}}^{\infty}$ is a sequence
satisfying (2.26), (2.30), (2.32) and (2.33), where $m_{1}\geq 1$. Then,

$\displaystyle\left\\{\left.||(x^{N}(t))^{(l)}||_{\infty}\right|N\geq
N_{1},\,l=0,1,\cdots,m_{1}\right\\}$

is bounded. If $f(x)$ and $g(x)$ are $C^{m_{1}-1}$, then

$\displaystyle\left\\{\left.||(u^{N}(t))^{(l)}||_{\infty}\right|N\geq
N_{1},\,l=0,1,\cdots,m_{1}-1\right\\}$

is bounded.

Proof.  Consider $(x_{r}^{N}(t))^{(m_{1})}$. From Lemma 5, it is a polynomial
of degree less than or equal to $N-r-m_{1}+1$. Therefore,

$\displaystyle(x_{r}^{N}(t))^{(m_{1})}=\displaystyle\sum_{n=0}^{N-r-
m_{1}+1}a_{n}^{N}(m_{1})L_{n}(t)$

where $L_{n}(t)$ is the Legendre polynomial of degree $n$. It is known that
$|L_{n}(t)|\leq 1$. Therefore, (2.33) implies that
$||(x_{r}^{N}(t))^{(m_{1})}||_{\infty}$ is bounded by ${\boldsymbol{d}}$ for
all $N\geq N_{1}$. Now, let us consider $(x_{r}^{N}(t))^{(m_{1}-1)}$. From
(2.14) we have,

$\displaystyle(x_{r}^{N}(t))^{(m_{1}-1)}$ $\displaystyle=$
$\displaystyle(x_{r}^{N}(t))^{(m_{1}-1)}|_{t=-1}+\int_{0}^{t}(x_{r}^{N}(s))^{(m_{1})}ds$
$\displaystyle=$
$\displaystyle\left[\begin{array}[]{cccccccccccccc}1&0&\cdots&0\end{array}\right]D^{m_{1}-1}(\bar{x}_{r}^{N})^{T}+\int_{0}^{t}(x_{r}^{N}(s))^{(m_{1})}ds$

So, $||(x_{r}^{N}(t))^{(m_{1}-1)}||_{\infty}$, $N\geq N_{1}$, is bounded
because of (2.32). Similarly, we can prove all derivatives of $x_{r}^{N}(t)$
of order less than $m_{1}$ are bounded. The same approach can also be applied
to prove the bound

$u^{N}(t)=\displaystyle\frac{\dot{x}_{r}^{N}(t)-f(x^{N}(t))}{g(x^{N}(t))}$

Because $f(x)$ and $g(x)$ have continuous derivatives of order less than or
equal to $m_{1}-1$, the boundedness of

$\left\\{\left.||(u^{N}(t))^{(l)}||_{\infty}\right|N\geq
N_{1},\,j=0,1,\cdots,m_{1}-1\right\\}$

follows the boundedness of $(x_{r}^{N}(t))^{(l)}$ proved above. $\Box$

Given any function $h(t)$ defined on $[-1,1]$. In the following, $U(h)$
represents an upper bound of $h(t)$ and $V(h)$ represents the total variation.

###### Lemma 7

Let $(x(t),u(t))$ be a solution of the differential equation (2.6). Suppose
$x_{r}^{(m)}(t)$ has bounded variation for some $m\geq 2$. Let $m_{1}$ be an
integer satisfying $1\leq m_{1}\leq m-1$. Then, there exist constants $M>0$
and $N_{1}>0$ so that for each integer $N\geq N_{1}$ the differential equation
(2.6) has a solution $(x^{N}(t),u^{N}(t))$ in which $x^{N}(t)$ consists of
polynomials of degree less than or equal to $N$. Furthermore, the pair
$(x^{N}(t),u^{N}(t))$ satisfies

$\displaystyle||x_{i}^{N}(t)-x_{i}(t)||_{\infty}$ $\displaystyle\leq$
$\displaystyle\displaystyle\frac{M||x_{r}||_{W^{m,2}}}{(N-r-
m_{1}+1)^{(m-m_{1})-3/4}},\;\;\;i=1,2,\cdots,r$ (3.23)
$\displaystyle||(x_{r}^{N}(t))^{(l)}-(x_{r}(t))^{(l)}||_{\infty}$
$\displaystyle\leq$
$\displaystyle\displaystyle\frac{M||x_{r}||_{W^{m,2}}}{(N-r-
m_{1}+1)^{(m-m_{1})-3/4}},\;\;\;\;l=1,2,\cdots,m_{1}$ (3.24)
$\displaystyle||u^{N}(t)-u(t)||_{\infty}$ $\displaystyle\leq$
$\displaystyle\displaystyle\frac{M||x_{r}||_{W^{m,2}}}{(N-r-
m_{1}+1)^{(m-m_{1})-3/4}}$ (3.25)

Furthermore, the spectral coefficients of $(x^{N}_{r})^{(m_{1})}(t)$ satisfy

$\displaystyle|a^{N}_{n}(m_{1})|\leq\displaystyle\frac{6(U(x^{(m_{1}+1)}_{r})+V(x^{(m_{1}+1)}_{r}))}{\sqrt{\pi}n^{3/2}},\;\;n=1,2,\cdots,N-r-1$
(3.26)

If $f(x)$ and $g(x)$ have Lipschitz continuous $L$th order partial derivatives
for some $L\leq m_{1}-1$, then

$\displaystyle||(u^{N}(t))^{(l)}-(u(t))^{(l)}||_{\infty}$ $\displaystyle\leq$
$\displaystyle\displaystyle\frac{M||x_{r}||_{W^{m,2}}}{(N-r-
m_{1}+1)^{(m-m_{1})-3/4}},\;\;\;l=1,\cdots,L$ (3.27)

Furthermore,

$\displaystyle\begin{array}[]{lll}x^{N}(-1)=x(-1)\\\ u^{N}(-1)=u(-1),&\mbox{
If }m_{1}\geq 2\end{array}$ (3.30)

###### Remark 3.3

In this lemma, if $x_{r}(t)$ has the $m$-th order strong derivative and if
$x_{r}^{(m)}(t)$ has bounded variation for some $m\geq 2$, then the
inequalities (3.23), (3.24), and (3.25) are slightly tighter.

$\displaystyle||x_{i}^{N}(t)-x_{i}(t)||_{\infty}$ $\displaystyle\leq$
$\displaystyle\displaystyle\frac{M||x_{r}||_{W^{m,2}}}{(N-r-
m_{1}+1)^{(m-m_{1})-1/2}},\;\;\;i=1,2,\cdots,r$ (3.31)
$\displaystyle||(x_{r}^{N}(t))^{(l)}-(x_{r}(t))^{(l)}||_{\infty}$
$\displaystyle\leq$
$\displaystyle\displaystyle\frac{M||x_{r}||_{W^{m,2}}}{(N-r-
m_{1}+1)^{(m-m_{1})-1/2}},\;\;\;\;l=1,2,\cdots,m_{1}$ (3.32)
$\displaystyle||u^{N}(t)-u(t)||_{\infty}$ $\displaystyle\leq$
$\displaystyle\displaystyle\frac{M||x_{r}||_{W^{m,2}}}{(N-r-
m_{1}+1)^{(m-m_{1})-1/2}}$ (3.33)

The proof is identical as that of Lemma 7 except that the error estimation in
(3) of Lemma 4 is used.

Proof.  Consider the Legendre series

$(x_{r})^{(m_{1})}(t)\sim\displaystyle\sum_{n=0}^{N-r-
m_{1}+1}a_{n}^{N}(m_{1})L_{n}(t)$

A sequence of polynomials $x_{1}^{N}(t),\cdots,x_{r+m_{1}}^{N}(t)$ is defined
as follows,

$\displaystyle x_{r+m_{1}}^{N}(t)$ $\displaystyle=$
$\displaystyle\displaystyle\sum_{n=0}^{N-r-m_{1}+1}a_{n}^{N}(m_{1})L_{n}(t)$
$\displaystyle x_{r+m_{1}-1}^{N}(t)$ $\displaystyle=$
$\displaystyle(x_{r})^{(m_{1}-1)}(-1)+\displaystyle{\int_{-1}^{t}}x_{r+m_{1}}^{N}(s)ds$
$\displaystyle\vdots$ $\displaystyle x_{r+1}^{N}(t)$ $\displaystyle=$
$\displaystyle\dot{x}_{r}(-1)+\displaystyle{\int_{-1}^{t}}x_{r+2}^{N}(s)ds$

and

$\displaystyle x_{i}^{N}(t)$ $\displaystyle=$ $\displaystyle
x_{i}(-1)+\displaystyle{\int_{-1}^{t}}x_{i+1}^{N}(s)ds,\;\;\mbox{ for }1\leq
i\leq r$

Define

$x^{N}(t)=\left[\begin{array}[]{cccccccccccccc}x_{1}^{N}(t)&\cdots&x^{N}_{r}(t)\end{array}\right]^{T}$

and define

$u^{N}(t)=\displaystyle\frac{x_{r+1}^{N}(t)-f(x^{N}(t))}{g(x^{N}(t))}$

From the definition of $x^{N}(t)$, we have $x^{N}(-1)=x(-1)$. If $m_{1}\geq
2$, then $x_{r+1}(-1)=\dot{x}_{r}(-1)$. From the definition of $u^{N}(t)$, we
know $u^{N}(-1)=u(-1)$ provided $m_{1}\geq 2$. Therefore,
$(x^{N}(t),u^{N}(t))$ satisfies (3.30). It is obvious that $x_{i}^{N}(t)$ is a
polynomial of degree less than or equal to $N$; and $(x^{N}(t),u^{N}(t))$
satisfies the differential equation (2.6). Because we assume
$V(x_{r}^{(m)})<\infty$, we have $x_{r}^{(m)}\in L^{2}$. From Lemma 4

$\displaystyle||x_{r+m_{1}}^{N}(t)-x_{r}^{(m_{1})}(t)||_{\infty}$
$\displaystyle=$
$\displaystyle||x_{r}^{(m_{1})}(t)-\displaystyle\sum_{n=0}^{N-r-
m_{1}+1}a_{n}^{N}(m_{1})L_{n}(t)||_{\infty}$ $\displaystyle\leq$
$\displaystyle C_{1}||x_{r}||_{W^{m,2}}(N-r-m_{1}+1)^{-(m-m_{1})+3/4}$

for some constant $C_{1}>0$. Therefore,

$\displaystyle|x^{N}_{r+m_{1}-1}(t)-(x_{r})^{(m_{1}-1)}(t)|$
$\displaystyle\leq$
$\displaystyle\displaystyle{\int_{-1}^{t}}|x^{N}_{r+m_{1}}(s)-(x_{r})^{(m_{1})}(s)|ds$
$\displaystyle\leq$ $\displaystyle 2C_{1}||x_{r}||_{W^{m,2}}(N-r-
m_{1}+1)^{-(m-m_{1})+3/4}$

Similarly, we can prove (3.23) and (3.24).

To prove (3.26), note that the spectral coefficient $a^{N}_{n}(m_{1})$ of
$(x_{r}^{N})^{(m_{1})}(t)$ is the same as the spectral coefficients of
$(x_{r})^{(m_{1})}(t)$. From Jackson’s Theorem (Lemma 2), we have

$|a^{N}_{n}(m_{1})|<\displaystyle\frac{6}{\sqrt{\pi}}(U(x^{(m_{1}+1)}_{r})+V(x^{(m_{1}+1)}_{r}))\displaystyle\frac{1}{n^{3/2}}$

In a bounded set around $x(t)$, we have $g(x)>\alpha>0$ for some $\alpha>0$
because $f$ and $g$ are Lipschitz continuous (Definition of Problem B).
Therefore, the function

$\displaystyle\frac{s-f(x)}{g(x)}$

is Lipschitz in a neighborhood of $(x,s)$, i.e. there exists a constant
$C_{2}$ independent of $N$ such that

$\displaystyle|u^{N}(t)-u(t)|$ $\displaystyle=$
$\displaystyle\left|\frac{x_{r+1}^{N}(t)-f(x^{N}(t))}{g(x^{N}(t))}-\frac{\dot{x}_{r}(t)-f(x(t))}{g(x(t))}\right|$
(3.34) $\displaystyle\leq$ $\displaystyle
C_{2}(|x_{r+1}^{N}(t)-\dot{x}_{r}(t)|+|x_{1}^{N}(t)-x_{1}(t)|+\cdots+|x_{r}^{N}(t)-x_{r}(t)|)$

Hence, (3.25) follows (3.23), (3.24) and (3.34) when $l=0$. Similarly, we can
prove (3.27) for $l\leq L$. $\Box$

Now, only after this lengthy work of preparation, we are ready to prove
Theorem 1.

Proof of Theorem 1: Let $(x^{\ast}(t),u^{\ast}(t))$ be an optimal solution to
Problem B. According to Lemma 7 and Remark 3.3, for any positive integer $N$
that is large enough, there exists a pair of functions
$(\hat{x}^{N}(t),\hat{u}^{N}(t))$ in which $\hat{x}^{N}(t)$ consists of
polynomials of degree less than or equal to $N$. Furthermore, the pair
satisfies the differential equation with initial conditions in Problem B and

$\displaystyle||\hat{x}^{N}(t)-x^{\ast}(t)||_{\infty}$ $\displaystyle<$
$\displaystyle\displaystyle\frac{L}{(N-r-m_{1}+1)^{m-m_{1}-1/2}}$ (3.35)
$\displaystyle||\hat{u}^{N}(t)-u^{\ast}(t)||_{\infty}$ $\displaystyle<$
$\displaystyle\displaystyle\frac{L}{(N-r-m_{1}+1)^{m-m_{1}-1/2}}$ (3.36)
$\displaystyle||(\hat{x}_{r}^{N}(t))^{(l)}-(x^{\ast}_{r}(t))^{(l)}||_{\infty}$
$\displaystyle<$ $\displaystyle\displaystyle\frac{L}{(N-r-
m_{1}+1)^{m-m_{1}-1/2}},\;\;\;1\leq l\leq m_{1}$ (3.37)

If we define

$\hat{\bar{u}}^{Nk}=\hat{u}^{N}(t_{k}),\hat{\bar{x}}^{Nk}=\hat{x}^{N}(t_{k})$

Then $\\{(\hat{\bar{x}}^{N},\hat{\bar{u}}^{N})\\}$ satisfies (2.26) and (2.27)
(Lemma 5 and 7). Because $\hat{x}_{r}^{N}(t)$ is a polynomial of degree less
than or equal to $N$ and because of (2.14), we know
$(\hat{x}^{N}_{r}(t))^{(j)}$ equals the interpolation polynomial of
$\hat{\bar{x}}^{N}_{r}(D^{T})^{j}$. So,

$\left[\begin{array}[]{cccccccccccccc}1&0&\cdots&0\end{array}\right]D^{j}(\hat{\bar{x}}_{r}^{N})^{T}=(\hat{x}^{N}_{r}(t))^{(j)}|_{t=-1}$

Therefore, (3.37) implies (2.32) if the bounds
$\underline{{\boldsymbol{b}}}_{j}$ and ${\boldsymbol{b}}_{j}$ are large
enough. In addition, the spectral coefficients of
$\hat{\bar{x}}^{N}_{r}(D^{T})^{m_{1}}$ is the same as the spectral
coefficients of $(\hat{x}^{N}_{r}(t))^{(m_{1})}$. From (3.26), (3.7) and
(3.6), we have

$\displaystyle\sum_{n=0}^{N-r-m_{1}+1}|a^{N}_{n}(m_{1})|\leq{\boldsymbol{d}}$

So, the spectral coefficients of $(\hat{x}_{r}^{N})^{(m_{1})}$ satisfies
(2.33). Because we select $\underline{{\boldsymbol{b}}}$ and
$\bar{\boldsymbol{b}}$ large enough so that the optimal trajectory of the
original continuous-time problem is contained in the interior of the region,
then (3.35) and (3.36) imply (2.30) for $N$ large enough. In summary, we have
proved that $(\hat{\bar{x}}^{N},\hat{\bar{u}}^{N})$ is a discrete feasible
trajectory satisfying all constraints, (2.26)-(2.33), in Problem ${\rm B}^{\rm
N}$.

Given any bounded control input $u(\cdot)$, because the system is globally
Lipschitz, it uniquely determines the trajectory $x(\cdot)$ if the initial
state is fixed. Therefore, the cost $J(x^{\ast}(\cdot),u^{\ast}(\cdot))$ can
be considered as a functional, denoted by ${\cal J}(u)$. Because all the
functions in Problem B are $C^{m}$ with $m\geq 2$, we know that ${\cal J}(u)$
has second order Fréchet derivative. By Lemma 1

$\displaystyle|J(x^{\ast}(\cdot),u^{\ast}(\cdot))-J(\hat{x}^{N}(\cdot),\hat{u}^{N}(\cdot))|$
(3.38) $\displaystyle=$ $\displaystyle|{\cal J}(u^{\ast})-{\cal
J}(\hat{u}^{N})|$ $\displaystyle\leq$ $\displaystyle
C_{1}(||u^{\ast}-\hat{u}^{N}||_{W^{m_{1}-1,\infty}}^{2})$ $\displaystyle\leq$
$\displaystyle\displaystyle\frac{C_{2}}{(N-r-m_{1}+1)^{2m-2m_{1}-1}}$

for some constant numbers $C_{1}$ and $C_{2}$. The last inequality is from
(3.36).

Now, consider $F(\hat{x}^{N}(t),\hat{u}^{N}(t))$ as a function of $t$. Let
$F^{N}(t)$ represent the polynomial interpolation of this function at
$t=t_{0},t_{1},\cdots,t_{N}$. Let $\hat{p}(t)$ be the best polynomial
approximation of $F(\hat{x}^{N}(t),\hat{u}^{N}(t))$ under the norm of
$L^{\infty}[-1,1]$. Then we have

$\displaystyle|J(\hat{x}^{N}(\cdot),\hat{u}^{N}(\cdot))-\bar{J}^{N}(\hat{\bar{x}}^{N},\hat{\bar{u}}^{N})|$
(3.39) $\displaystyle=$
$\displaystyle|J(\hat{x}^{N}(\cdot),\hat{u}^{N}(\cdot))-\displaystyle\sum_{k=0}^{N}F(\hat{\bar{x}}^{Nk},\hat{\bar{u}}^{Nk})w_{k}-E(\hat{\bar{x}}^{N0},\hat{\bar{x}}^{NN})|$
$\displaystyle=$
$\displaystyle\left|\int_{-1}^{1}F(\hat{x}^{N}(t),\hat{u}^{N}(t))dt-\int_{-1}^{1}F^{N}(t)dt\right|$
$\displaystyle\leq$
$\displaystyle\int_{-1}^{1}|F(\hat{x}^{N}(t),\hat{u}^{N}(t))-F^{N}(t)|dt$
$\displaystyle\leq$ $\displaystyle
2(1+\Lambda_{N})||\hat{p}(t)-F(\hat{x}^{N}(t),\hat{u}^{N}(t))||_{\infty}$

where

$\begin{array}[]{llllllllll}\Lambda_{N}\leq\displaystyle\frac{2}{\pi}log(N+1)+0.685\cdots\end{array}$
(3.40)

is the Lebesgue constant. The inequality (3.39) is a corollary of Lemma 3.
Because $f(\cdot)$, $g(\cdot)$, and $F(\cdot)$ are $C^{m}$, it is known (Lemma
4) that the best polynomial approximation satisfies

$||\hat{p}(t)-F(\hat{x}^{N}(t),\hat{u}^{N}(t))||_{\infty}\leq\displaystyle\frac{C_{3}}{N^{m_{1}-1}}||F(\hat{x}^{N}(t),\hat{u}^{N}(t))||_{W^{m_{1}-1,\infty}}$

Because of Lemma 6,
$\\{||F(\hat{x}^{N}(t),\hat{u}^{N}(t))||_{W^{m_{1}-1,\infty}}|N\geq N_{1}\\}$
is bounded. Therefore,

$\displaystyle|J(\hat{x}^{N}(\cdot),\hat{u}^{N}(\cdot))-\bar{J}^{N}(\hat{\bar{x}}^{N},\hat{\bar{u}}^{N})|\leq\displaystyle\frac{(1+\Lambda_{N})C_{4}}{N^{m_{1}-1}}\leq\displaystyle\frac{C_{5}}{N^{\alpha}}$
(3.41)

for some constant numbers $C_{4}$ and $C_{5}$ independent of $N$ and any
$\alpha<m_{1}-1$. Let

$\begin{array}[]{llllllllll}\\{(\bar{x}^{\ast N},\bar{u}^{\ast
N})\\}_{N=N_{0}}^{\infty}\end{array}$ (3.42)

be a sequence of optimal discrete solutions. Its interpolation is denoted by
$(x^{\ast N}(t),u^{\ast N}(t))$. Then, similar to the derivation above, we can
prove

$\displaystyle|J(x^{\ast N}(\cdot),u^{\ast
N}(\cdot))-\bar{J}^{N}(\bar{x}^{\ast N},\bar{u}^{\ast N})|$
$\displaystyle\leq$ $\displaystyle 2(1+\Lambda_{N})||p^{N}(t)-F(x^{\ast
N}(t),u^{\ast N}(t))||_{\infty}$ $\displaystyle\leq$
$\displaystyle\displaystyle\frac{C_{6}(1+\Lambda_{N})}{N^{m_{1}-1}}||F(x^{\ast
N}(t),u^{\ast N}(t))||_{W^{m_{1}-1,\infty}}$

where $p^{N}(t)$ is the best polynomial approximation of $F(x^{\ast
N}(t),u^{\ast N}(t))$ with degree less than or equal to $N$. Because of Lemma
6, $||F(x^{\ast N}(t),u^{\ast N}(t))||_{W^{m_{1}-1,\infty}}|N\geq N_{1}\\}$ is
bounded. So

$\displaystyle|J(x^{\ast N}(\cdot),u^{\ast
N}(\cdot))-\bar{J}^{N}(\bar{x}^{\ast N},\bar{u}^{\ast
N})|\leq\displaystyle\frac{C_{7}}{N^{\alpha}}$ (3.44)

for some constant $C_{7}>0$. Now, we are ready to piece together the puzzle of
inequalities and finalize the proof.

$\begin{array}[]{rcllll}&&J(x^{\ast}(\cdot),u^{\ast}(\cdot))\\\
&\leq&J(x^{\ast N}(\cdot),u^{\ast N}(\cdot))&\left(\begin{array}[]{ll}(x^{\ast
N}(t),u^{\ast N}(t))\mbox{ is a feasible }\\\ \mbox{trajectory (Lemma
\ref{lemma1})}\end{array}\right)\\\ &\leq&\bar{J}^{N}(\bar{x}^{\ast
N},\bar{u}^{\ast N})+\displaystyle\frac{C_{7}}{N^{\alpha}}&\left(\mbox{
inequality }(\ref{eq3_11})\right)\\\
&\leq&\bar{J}^{N}(\hat{x}^{N},\hat{u}^{N})+\displaystyle\frac{C_{7}}{N^{\alpha}}&\left(\begin{array}[]{ll}(\hat{x}^{N},\hat{u}^{N})\mbox{
is a feasible discrete }\\\ \mbox{trajectory and }(\bar{x}^{\ast
N},\bar{u}^{\ast N})\mbox{ is optimal}\end{array}\right)\\\
&\leq&J(\hat{x}^{N}(\cdot),\hat{u}^{N}(\cdot))+\displaystyle\frac{C_{5}}{N^{\alpha}}+\displaystyle\frac{C_{7}}{N^{\alpha}}&\left(\mbox{
inequality }(\ref{eq3_8})\right)\\\
&\leq&J(x^{\ast}(\cdot),u^{\ast}(\cdot))+\displaystyle\frac{C_{2}}{(N-r-
m_{1}-1)^{2m-2m_{1}-1}}\\\
&&+\displaystyle\frac{C_{5}}{N^{\alpha}}+\displaystyle\frac{C_{7}}{N^{\alpha}}&\left(\mbox{
inequality }(\ref{eq3_6})\right)\end{array}$

Therefore,

$0\leq J(x^{\ast N}(\cdot),u^{\ast
N}(\cdot))-J(x^{\ast}(\cdot),u^{\ast}(\cdot))\leq\displaystyle\frac{C_{2}}{(N-r-
m_{1}-1)^{2m-2m_{1}-1}}+\displaystyle\frac{C_{5}}{N^{\alpha}}+\displaystyle\frac{C_{7}}{N^{\alpha}}$

This inequality implies (3.9). Furthermore, (3.9) and (3.44) imply (3.10).
$\Box$

According to Theorem 1, the convergence rate of the approximate cost is
determined by two terms with the rates

$\displaystyle\displaystyle\frac{1}{(N-r-
m_{1}-1)^{2m-2m_{1}-1}}\sim\displaystyle\frac{1}{N^{2m-2m_{1}-1}}$ (3.45)

and

$\displaystyle\displaystyle\frac{1}{N^{\alpha}}\sim\displaystyle\frac{1}{N^{m_{1}-1}}$
(3.46)

where $m$, the smoothness of $x^{\ast}(t)$, is fixed. However, $m_{1}$ can be
selected provided $f(\cdot)$, $g(\cdot)$, and $F(\cdot)$ are smooth enough.
Note that increasing $m_{1}$ will increase the rate defined by (3.46), but
decrease the rate defined by (3.45). There is a value of $m_{1}$ that
determines the maximum rate. Given any real number $a\in\Re$, let $[a]$ be the
greatest integer less than or equal to $a$.

###### Corollary 1

Under the same assumption as Theorem 1, the convergence rate of $J(x^{\ast
N}(\cdot),u^{\ast N}(\cdot))$ and $\bar{J}^{N}(\bar{x}^{\ast N},\bar{u}^{\ast
N})$ is

$O\left(\frac{1}{N^{[\frac{2m}{3}]-\delta}}\right)$

in which

$\delta=\left\\{\begin{array}[]{ll}1&0<\gamma<\frac{2}{3}\\\
3\left(1-\gamma\right)&\gamma\geq\frac{2}{3}\\\ 1-\mbox{any positive
number},&\gamma=0\end{array}\right.$

where $\gamma=\frac{2m}{3}-[\frac{2m}{3}]$. To achieve this rate,

$m_{1}=\left\\{\begin{array}[]{lll}\left[\displaystyle\frac{2m}{3}\right],&0\leq\gamma<\frac{2}{3}\\\
\\\
\left[\displaystyle\frac{2m}{3}\right]+1,&\gamma\geq\frac{2}{3}\end{array}\right.$

Proof.  The optimal convergence rate is determined by

$\max_{2\leq m_{1}\leq m-1}\min\\{2m-2m_{1}-1,\,m_{1}-1\\}$

The maxmin is achieved at

$m_{1}=\frac{2m}{3}$

However, it may not be an integer. If $m_{1}$ is not an integer, we have two
options,

$m_{1}=[\frac{2m}{3}]\mbox{ or }[\frac{2m}{3}]+1$

If we define

$\gamma=\frac{2m}{3}-[\frac{2m}{3}]$

then either $m_{1}=\frac{2m}{3}-\gamma$ or $m_{1}=\frac{2m}{3}-\gamma+1$. It
is straightforward to verify that

$\min\\{2m-2m_{1}-1,\,m_{1}-1\\}=\left\\{\begin{array}[]{lll}\frac{2m}{3}-\gamma-1,&m_{1}=\frac{2m}{3}-\gamma\\\
\frac{2m}{3}-\gamma-3(1-\gamma),&m_{1}=\frac{2m}{3}-\gamma+1\end{array}\right.$

Therefore, $\frac{2m}{3}-1-\gamma$ is larger when $\gamma<\frac{2}{3}$, and
$\frac{2m}{3}-1-2(1-\gamma)$ is larger if $\gamma\geq\frac{2}{3}$. The special
case at $\gamma=0$ is because of (3.46) when $\frac{2m}{3}$ equals an integer.
$\Box$

Different from numerical computations of differential equations, solving an
optimal control problem requires the approximation, (2.1), of the integration
as an addition to the approximation, (2.26), of the differential equation. The
contributions of these approximations to the overall approximation error are
different; and the errors are inversely related to each other. The following
theorem indicates that the rate (3.45) is due to the approximation error of
the differential equation and the rate (3.46) is due to the approximation
error of the quadrature integration rule (4.9). To verify this fact, we define
the following discretization problem with exact integration.

Problem ${\rm B}^{\rm N}$(J) Find $\bar{x}^{Nk}\in\Re^{r}$ and
$\bar{u}^{Nk}\in\Re$, $k\ =\ 0,1,\ldots,N$, that minimize

$\displaystyle J(x^{N}(\cdot),u^{N}(\cdot))$ $\displaystyle=$
$\displaystyle\int_{-1}^{1}F(x^{N}(t),u^{N}(t))\ dt+E(x^{N}(-1),x^{N}(1))$
(3.47)

subject to

$\displaystyle\left\\{\begin{array}[]{rcl}D(\bar{x}_{1}^{N})^{T}&=&(\bar{x}_{2}^{N})^{T}\\\
D(\bar{x}_{2}^{N})^{T}&=&(\bar{x}_{3}^{N})^{T}\\\ &\vdots&\\\
D(\bar{x}_{r-1}^{N})^{T}&=&(\bar{x}_{r}^{N})^{T}\\\
D(\bar{x}_{r}^{N})^{T}&=&\left[\begin{array}[]{cccccccccccccc}f(\bar{x}^{N0})+g(\bar{x}^{N0})\bar{u}^{N0}\\\
\vdots\\\ f(\bar{x}^{NN})+g(\bar{x}^{NN})\bar{u}^{NN}\end{array}\right]\\\
\end{array}\right.$ (3.56) $\displaystyle\bar{x}^{N0}=x_{0}$ (3.57)
$\displaystyle\underline{{\boldsymbol{b}}}\leq\left[\begin{array}[]{cccccccccccccc}\bar{x}^{Nk}\\\
\bar{u}^{Nk}\end{array}\right]\ \leq\ \bar{\boldsymbol{b}},\;\;\;\;\mbox{ for
all }0\leq k\leq N$ (3.60)
$\displaystyle\underline{{\boldsymbol{b}}}_{j}\leq\left[\begin{array}[]{cccccccccccccc}1&0&\cdots&0\end{array}\right]D^{j}(\bar{x}_{r}^{N})^{T}\
\leq\bar{\boldsymbol{b}}_{j},\;\;1\leq j\leq m_{1}-1$ (3.62)
$\displaystyle\displaystyle\sum_{n=0}^{N-r-
m_{1}+1}|a^{N}_{n}(m_{1})|\leq{\boldsymbol{d}}$ (3.63)

In Problem ${\rm B}^{\rm N}$(J), $(x^{N}(t),u^{N}(t))$ is the interpolation of
$(\bar{x}^{N},\bar{u}^{N})$. In this discretization, we approximate the
differential equation by the PS method. However, the integration in the cost
function is exact. In this case, the overall error is controlled by the single
rate (3.45) rather than the two-rate convergence of Problem ${\rm B}^{\rm N}$.
Without the integration error of the cost function, the convergence rate is
improved to $\displaystyle\frac{1}{N^{2m-3}}$; and the smoothness requirement
can be reduced to $m\geq 2$.

###### Theorem 2

Suppose Problem B has an optimal solution $(x^{\ast}(t),u^{\ast}(t))$ in which
the strong derivative $(x_{r}^{\ast}(t))^{(m)}$ has bounded variation for some
$m\geq 2$. In Problem ${\rm B}^{\rm N}$(J), select $m_{1}$ so that $1\leq
m_{1}\leq m-1$. Suppose $f(\cdot)$, $g(\cdot)$, and $F(\cdot)$ are $C^{m}$.
Suppose all other bounds in Problem ${\rm B}^{\rm N}$ are large enough. Given
any sequence

$\begin{array}[]{llllllllll}\\{(\bar{x}^{\ast N},\bar{u}^{\ast N})\\}_{N\geq
N_{1}}\end{array}$ (3.64)

of optimal solutions of Problem ${\rm B}^{\rm N}$(J). Then the cost of (3.64)
converges to the optimal cost at the following rate

$\displaystyle\left|J(x^{\ast}(\cdot),u^{\ast}(\cdot))-J(x^{\ast
N}(\cdot),u^{\ast N}(\cdot))\right|\leq\displaystyle\frac{M_{1}}{(N-r-
m_{1}-1)^{2m-2m_{1}-1}}$ (3.65)

for some constants $M_{1}$ independent of $N$.

Proof.  Let

$\begin{array}[]{llllllllll}\\{(\bar{x}^{\ast N},\bar{u}^{\ast
N})\\}_{N=N_{0}}^{\infty}\end{array}$ (3.66)

be a sequence of optimal solutions of Problem ${\rm B}^{\rm N}$(J). According
to Lemma 7 and Remark 3.3, for any positive integer $N$ that is large enough,
there exists a pair of functions $(\hat{x}^{N}(t),\hat{u}^{N}(t))$ in which
$\hat{x}^{N}(t)$ consists of polynomials of degree less than or equal to $N$.
Furthermore, the pair satisfies the differential equation in Problem B and the
inequalities (3.35), (3.36), and (3.37). If we define

$\hat{\bar{u}}^{Nk}=\hat{u}^{N}(t_{k}),\hat{\bar{x}}^{Nk}=\hat{x}^{N}(t_{k})$

Then, from the first part in the proof of Theorem 1,
$(\hat{\bar{x}}^{N},\hat{\bar{u}}^{N})$ is a discrete feasible solution
satisfying all constraints in Problem ${\rm B}^{\rm N}$(J). By Lemma 1

$\displaystyle|J(x^{\ast}(\cdot),u^{\ast}(\cdot))-J(\hat{x}^{N}(\cdot),\hat{u}^{N}(\cdot))|$
(3.67) $\displaystyle=$ $\displaystyle|{\cal J}(u^{\ast})-{\cal
J}(\hat{u}^{N})|$ $\displaystyle\leq$ $\displaystyle
C_{1}(||u^{\ast}-\hat{u}^{N}||_{W^{m_{1}-1,\infty}}^{2})$ $\displaystyle\leq$
$\displaystyle\displaystyle\frac{C_{2}}{(N-r-m_{1}+1)^{2m-2m_{1}-1}}$

for some constant numbers $C_{1}$ and $C_{2}$. The last inequality is from
(3.25). The interpolation $(x^{\ast N}(t),u^{\ast N}(t))$ of (3.66) is a
feasible trajectory of Problem B (Lemma 5). Thus,

$\begin{array}[]{rcllll}&&J(x^{\ast}(\cdot),u^{\ast}(\cdot))\\\
&\leq&J(x^{\ast N}(\cdot),u^{\ast N}(\cdot))&\left(\begin{array}[]{ll}(x^{\ast
N}(t),u^{\ast N}(t))\mbox{ is a feasible }\\\ \mbox{trajectory (Lemma
\ref{lemma1})}\end{array}\right)\\\
&\leq&J(\hat{x}^{N}(\cdot),\hat{u}^{N}(\cdot))&\left(\begin{array}[]{ll}(\hat{\bar{x}}^{N},\hat{\bar{u}}^{N})\mbox{
is a feasible discrete }\\\ \mbox{trajectory and }(\bar{x}^{\ast
N},\bar{u}^{\ast N})\mbox{ is optimal}\end{array}\right)\\\
&\leq&J(x^{\ast}(\cdot),u^{\ast}(\cdot))+\displaystyle\frac{C_{2}}{(N-r-
m_{1}-1)^{2m-2m_{1}-1}}&\left(\mbox{ inequality
}(\ref{eq3_6a})\right)\end{array}$

Therefore,

$0\leq J(x^{\ast N}(\cdot),u^{\ast
N}(\cdot))-J(x^{\ast}(\cdot),u^{\ast}(\cdot))\leq\displaystyle\frac{C_{2}}{(N-r-
m_{1}-1)^{2m-2m_{1}-1}}$

$\Box$

##  4 Existence and Convergence of Approximate Optimal Solutions

In Section 3, the rate of convergence for the cost function is proved.
However, the results do not guarantee the convergence of the approximate
optimal trajectory $\\{(x^{N}(t),u^{N}(t))\\}$. In this section, we prove the
existence of feasible trajectories for Problem ${\rm B}^{\rm N}$ and the
existence of a convergent subsequence in any set of approximate optimal
solutions. In addition, we consider a larger family of problems. Different
from Section 2 where Problem B does not contain constraints other than the
control system, in this section the problem of optimal control may contain
nonlinear path constraints. Furthermore, general endpoint conditions are
allowed, rather than being limited to the initial value problem as in the
previous sections.

Problem B: Determine the state-control function pair $(x(t),u(t))$,
$x\in\Re^{r}$ and $u\in\Re$, that minimizes the cost function

$\displaystyle J(x(\cdot),u(\cdot))$ $\displaystyle=$
$\displaystyle\int_{-1}^{1}F(x(t),u(t))\ dt+E(x(-1),x(1))$ (4.1)

subject to the state equation

$\displaystyle\left\\{\begin{array}[]{lll}\dot{x}_{1}=x_{2}\\\ \;\;\;\vdots\\\
\dot{x}_{r-1}=x_{r}\\\ \dot{x}_{r}=f(x)+g(x)u\end{array}\right.$ (4.6)

end-point conditions

$\displaystyle e(x(-1),x(1))$ $\displaystyle=$ $\displaystyle 0$ (4.7)

and state-control constraints

$\displaystyle h(x(t),u(t))$ $\displaystyle\leq$ $\displaystyle 0$ (4.8)

where $x\in\Re^{r}$, $u\in\Re$, and $F:\Re^{r}\times\Re\to\Re$,
$E:\Re^{r}\times\Re^{r}\to\Re$, $f:\Re^{r}\to\Re$, $g:\Re^{r}\to\Re$
$e:\Re^{r}\times\Re^{r}\to\Re^{N_{e}}$ and $h:\Re^{r}\times\Re\to\Re^{N_{h}}$
are all Lipschitz continuous functions with respect to their arguments. In
addition, we assume $g(x)\neq 0$ for all $x$. The corresponding discretization
is defined as follows.

Problem ${\bf B}^{\bf N}$: Find $\bar{x}^{Nk}\in\Re^{r}$ and
$\bar{u}^{Nk}\in\Re$, $k\ =\ 0,1,\ldots,N$, that minimize

$\displaystyle\bar{J}^{N}(\bar{x}^{N},\bar{u}^{N})$ $\displaystyle=$
$\displaystyle\sum_{k=0}^{N}F(\bar{x}^{Nk},\bar{u}^{Nk})w_{k}+E(\bar{x}^{N0},\bar{x}^{NN})$
(4.9)

subject to

$\displaystyle\left\\{\begin{array}[]{rcl}D(\bar{x}_{1}^{N})^{T}&=&(\bar{x}_{2}^{N})^{T}\\\
D(\bar{x}_{2}^{N})^{T}&=&(\bar{x}_{3}^{N})^{T}\\\ &\vdots&\\\
D(\bar{x}_{r-1}^{N})^{T}&=&(\bar{x}_{r}^{N})^{T}\\\
D(\bar{x}_{r}^{N})^{T}&=&\left[\begin{array}[]{cccccccccccccc}f(\bar{x}^{N0})+g(\bar{x}^{N0})\bar{u}^{N0}\\\
\vdots\\\ f(\bar{x}^{NN})+g(\bar{x}^{NN})\bar{u}^{NN}\end{array}\right]\\\
\end{array}\right.$ (4.18)
$\displaystyle\|e(\bar{x}^{N0},\bar{x}^{NN})\|_{\infty}\leq(N-r-1)^{-\beta}$
(4.19) $\displaystyle
h(\bar{x}^{Nk},\bar{u}^{Nk})\leq(N-r-1)^{-\beta}\cdot\mathbf{1},\qquad\ \ \ \
\mbox{ for all }0\leq k\leq N$ (4.20)
$\displaystyle\underline{{\boldsymbol{b}}}\leq\left[\begin{array}[]{cccccccccccccc}\bar{x}^{Nk}\\\
\bar{u}^{Nk}\end{array}\right]\ \leq\ \bar{\boldsymbol{b}},\;\;\;\;\mbox{ for
all }0\leq k\leq N$ (4.23)
$\displaystyle\underline{{\boldsymbol{b}}}_{j}\leq\left[\begin{array}[]{cccccccccccccc}1&0&\cdots&0\end{array}\right]D^{j}(\bar{x}_{r}^{N})^{T}\
\leq\bar{\boldsymbol{b}}_{j},\mbox{ if }1\leq j\leq m_{1}-1\mbox{ and
}m_{1}\geq 2$ (4.25) $\displaystyle\displaystyle\sum_{n=0}^{N-r-
m_{1}+1}|a^{N}_{n}(m_{1})|\leq{\boldsymbol{d}}$ (4.26)

The discretization is almost identical to the one used in the previous
sections except for the path constraints and the endpoint conditions, which
must be treated with care. Note that in Problem ${\rm B}^{\rm N}$the right
sides of (4.19) and (4.20) are not zero. It is necessary to relax (4.7) and
(4.8) by a small margin for the reason of feasibility. The margin approaches
zero as $N$ is increased. Without this relaxation, it is shown by a counter
example in [7] that Problem ${\rm B}^{\rm N}$ may have no feasible
trajectories.

Some feasibility and convergence results were proved in [7], which take the
form of consistent approximation theory based on the convergence assumption
about $\\{\dot{x}_{r}^{N}(t)\\}$ and $\\{\bar{x}^{N0}\\}$. The goal of this
section is to remove this bothersome assumption by using a fundamentally
different approach. In addition, the proofs in this section are not based on
necessary conditions of optimal control and any coercivity assumption, which
are widely used in existing work on the convergence of direct optimal control
methods. Before we introduce main results in this section, some useful results
from [7] are summarized in the following Lemma.

###### Lemma 8

([7]) Suppose Problem B has an optimal solution $(x^{\ast}(t),u^{\ast}(t))$
satisfying $x^{\ast}_{r}(t)\in W^{m,\infty}$, $m\geq 2$. Let
$\\{(\bar{x}^{N},\bar{u}^{N})\\}_{N=N_{1}}^{\infty}$ be a sequence of feasible
solutions to (4.18)-(4.23). Suppose there is a subsequence
$\\{N_{j}\\}_{j=1}^{\infty}$ of $\\{N\\}_{N=1}^{\infty}$ such that the
sequence $\left\\{\bar{x}^{N_{j}0}\right\\}_{j=1}^{\infty}$ converges as
$N_{j}\rightarrow\infty$. Suppose there exists a continuous function $q(t)$
such that $\dot{x}_{r}^{N_{j}}(t)$ converges to $q(t)$ uniformly in $[-1,1]$.
Then, there exists $(x^{\infty}(t),u^{\infty}(t))$ satisfying (4.6)-(4.8) such
that the following limits converge uniformly in $[-1,1]$.

$\displaystyle\lim_{N_{j}\rightarrow\infty}(x^{N_{j}}(t)-x^{\infty}(t))=0$
(4.27)
$\displaystyle\lim_{N_{j}\rightarrow\infty}(u^{N_{j}}(t)-u^{\infty}(t))=0$
(4.28)
$\displaystyle\lim_{N_{j}\rightarrow\infty}\bar{J}^{N_{j}}(\bar{x}^{N_{j}},\bar{u}^{N_{j}})=J(x(\cdot),u(\cdot))$
(4.29)
$\displaystyle\lim_{N_{j}\rightarrow\infty}J(x^{N_{j}},u^{N_{j}})=J(x(\cdot),u(\cdot))$
(4.30)

In addition to the above assumptions, if
$\\{(\bar{x}^{N},\bar{u}^{N})\\}_{N=N_{1}}^{\infty}$ is a sequence of optimal
solutions subject to the constraints (4.18)-(4.23), then
$(x^{\infty}(t),u^{\infty}(t))$ must be an optimal solution to Problem B.

The following are the two main theorems of this section. Relative to [14],
these results has a tightened bounds for $m_{1}$ and $\beta$.

###### Theorem 3

(Existence of solutions) Consider Problem B and Problem ${\rm B}^{\rm N}$
defined in Section 2. Suppose Problem B has a feasible trajectory
$(x(t),u(t))$ in which $(x_{r}(t))^{(m)}$ has bounded variation for some
$m\geq 2$. In Problem ${\rm B}^{\rm N}$, let $m_{1}$ be any integer and
$\beta$ be any real number satisfying $1\leq m_{1}\leq m-1$ and
$0<\beta<(m-m_{1})-\frac{3}{4}$. Then, there exists $N_{1}>0$ so that, for all
$N\geq N_{1}$, Problem ${\rm B}^{\rm N}$ has a feasible trajectory satisfying
(4.18)-(4.26). Furthermore, $(x(t),u(t))$ and the interpolation
$(x^{N}(t),u^{N}(t))$ satisfy (3.23)-(3.26).

###### Theorem 4

(Convergence) Consider Problem B and Problem ${\rm B}^{\rm N}$ defined in
Section 2. Suppose Problem B has an optimal solution
$(x^{\ast}(t),u^{\ast}(t))$ in which $(x^{\ast}_{r}(t))^{(m)}$ has bounded
variation for some $m\geq 3$. In Problem ${\rm B}^{\rm N}$, let $m_{1}$ be any
integer and $\beta$ be any real number satisfying $2\leq m_{1}\leq m-1$ and
$0<\beta<(m-m_{1})-\frac{3}{4}$. Then for any sequence $\\{(\bar{x}^{\ast
N},\bar{u}^{\ast N})\\}_{N=N_{1}}^{\infty}$ of optimal solutions of Problem
${\rm B}^{\rm N}$, there exists a subsequence, $\\{(\bar{x}^{\ast
N_{j}},\bar{u}^{\ast N_{j}})\\}_{j\geq 1}^{\infty}$, and an optimal solution,
$(x^{\ast}(t),u^{\ast}(t))$, of Problem B so that the following limits
converge uniformly in $[-1,1]$

$\displaystyle\lim_{N_{j}\rightarrow\infty}(x^{\ast N_{j}}(t)-x^{\ast}(t))$
$\displaystyle=$ $\displaystyle 0$
$\displaystyle\lim_{N_{j}\rightarrow\infty}(u^{\ast N_{j}}(t)-u^{\ast}(t))$
$\displaystyle=$ $\displaystyle 0$ (4.31)
$\displaystyle\lim_{N_{j}\rightarrow\infty}\bar{J}^{N_{j}}(\bar{x}^{\ast
N_{j}},\bar{u}^{\ast N_{j}})$ $\displaystyle=$ $\displaystyle
J(x^{\ast}(\cdot),u^{\ast}(\cdot))$
$\displaystyle\lim_{N_{j}\rightarrow\infty}J(x^{\ast N_{j}}(\cdot),u^{\ast
N_{j}}(\cdot))$ $\displaystyle=$ $\displaystyle
J(x^{\ast}(\cdot),u^{\ast}(\cdot))$

where $(x^{\ast N_{j}}(t),u^{\ast N_{j}}(t))$ is the interpolation of
$(\bar{x}^{\ast N},\bar{u}^{\ast N})$.

###### Remark 4.1

The integers $m$ and $m_{1}$ in Theorem 3 are smaller than those in Theorem 4,
i.e. the existence theorem is proved under a weaker smoothness assumption than
the convergence theorem.

To prove these theorems, we first briefly review some results on real analysis
and then prove a lemma. Given a sequence of functions
$\\{f_{k}(t)\\}_{k=1}^{\infty}$ defined on $[a,b]$. It is said to be uniformly
equicontinuous if for every $\epsilon>0$, there exists a $\delta>0$ such that
for all $t$, $t^{\prime}$ in $[a,b]$ with $|t^{\prime}-t|<\delta$, we have

$|f_{k}(t)-f_{k}(t^{\prime})|<\epsilon$

for all $k\geq 1$. The following Proposition and Theorem are standard in real
analysis [20].

###### Proposition 1

If $f_{k}(t)$ is differentiable for all $k$, and if
$\\{\dot{f}_{k}(t)\\}_{k=1}^{\infty}$ is bounded. Then,
$\\{f_{k}(t)\\}_{k=1}^{\infty}$ is uniformly equicontinuous.

###### Theorem 5

(Arzelà-Ascoli Theorem) Consider a sequence of continuous functions
$\\{h_{n}(t)\\}_{n=1}^{\infty}$ defined on a closed interval $[a,b]$ of the
real line with real values. If this sequence is uniformly bounded and
uniformly equicontinuous, then it admits a subsequence which converges
uniformly.

###### Lemma 9

([14])Let $\\{(\bar{x}^{N},\bar{u}^{N})\\}_{N=N_{1}}^{\infty}$ be a sequence
satisfying (4.18)-(4.23). Assume the set

$\displaystyle\left\\{\left.||\ddot{x}^{N}_{r}(t)||_{\infty}\right|N\geq
N_{1}\right\\}$ (4.32)

is bounded. Then, there exists $(x^{\infty}(t),u^{\infty}(t))$ satisfying
(4.6)-(4.8) and a subsequence
$\\{(\bar{x}^{N_{j}},\bar{u}^{N_{j}})\\}_{N_{j}\geq N_{1}}^{\infty}$ such that
(4.27), (4.28), (4.29) and (4.30) hold. Furthermore, if
$\\{(\bar{x}^{N},\bar{u}^{N})\\}_{N=N_{1}}^{\infty}$ is a sequence of optimal
solutions to Problem ${\rm B}^{\rm N}$, then $(x^{\infty}(t),u^{\infty}(t))$
must be an optimal solution to Problem B.

Proof.  Let $x_{r}^{N}(t)$ be the interpolation polynomial of
$\bar{x}_{r}^{N}$. Because (4.32) is a bounded set, we know that the sequence
of functions $\\{\dot{x}_{r}^{N}(t)|N\geq N_{1}\\}$ is uniformly
equicontinuous (Proposition 1). By the Arzelà-Ascoli Theorem, a subsequence
$\\{\dot{x}_{r}^{N_{j}}(t)\\}$ converges uniformly to a continuous function
$q(t)$. In addition, because of (4.23), we can select the subsequence so that
$\\{\bar{x}^{N_{j}0}\\}_{N_{j}\geq N_{1}}^{\infty}$ is convergent. Therefore,
all conclusions in Lemma 8 hold true. $\Box$

Now, we are ready to prove the theorems.

Proof of Theorem 3: For the feasible trajectory $(x(t),u(t))$, consider the
pair $(x^{N}(t),u^{N}(t))$ in Lemma 7 that satisfies the differential equation
(4.6). Define

$\displaystyle\begin{array}[]{rcl}\bar{x}^{Nk}&=&x^{N}(t_{k})\\\
\bar{u}^{Nk}&=&u^{N}(t_{k})\end{array}$ (4.35)

for $0\leq k\leq N$. From Lemma 5, we know that
$\\{(\bar{x}^{N},\bar{u}^{N})\\}$ satisfies the discrete equations in (4.18).
In the next we prove that the mixed state-control constraint (4.20) is
satisfied. Because $h$ is Lipschitz continuous and because of (3.23) and
(3.25), there exists a constant $C$ independent of $N$ so that

$\displaystyle\|h(x(t),u(t))-h(x^{N}(t),u^{N}(t))\|$ $\displaystyle\leq$
$\displaystyle
C(|x_{1}(t)-x^{N}_{1}(t)|+\cdots+|x_{r}(t)-x^{N}_{r}(t)|+|u(t)-u^{N}(t)|)$
$\displaystyle\leq$ $\displaystyle CMV||x_{r}||_{W^{m,2}}(r+1)(N-r-
m_{1}+1)^{-(m-m_{1})+3/4}$

Hence

$\displaystyle h(x^{N}(t),u^{N}(t))$ $\displaystyle\leq$ $\displaystyle
h(x(t),u(t))+CM||x_{r}||_{W^{m,2}}(r+1)(N-r-
m_{1}+1)^{-(m-m_{1})+3/4}\cdot\mathbf{1}$ $\displaystyle\leq$ $\displaystyle
CM||x_{r}||_{W^{m,2}}(r+1)(N-r-m_{1}+1)^{-(m-m_{1})+3/4}$

Because $\beta<m-m_{1}-\frac{3}{4}$, there exists a positive integer $N_{1}$
such that, for all $N>N_{1}$,

$\displaystyle CM||x_{r}||_{W^{m,2}}(r+1)(N-r-m_{1}+1)^{-(m-m_{1})+3/4}$
$\displaystyle\leq$ $\displaystyle(N-r-1)^{-\beta}$

Therefore $x^{N}_{1}(t_{k})$, $\ldots$, $x^{N}_{r}(t_{k})$, $u^{N}(t_{k})$,
$k=0,1,\ldots,N$, satisfy the mixed state and control constraint (4.20) for
all $N>N_{1}$.

By a similar procedure, we can prove that the endpoint condition (4.19) is
satisfied. Because $x_{r}^{N}(t)$ is a polynomial of degree less than or equal
to $N$, and because of (2.14) and (4.35), we know $(x^{N}_{r}(t))^{(j)}$
equals the interpolation polynomial of $\bar{x}^{N}_{r}(D^{T})^{j}$. So,

$\left[\begin{array}[]{cccccccccccccc}1&0&\cdots&0\end{array}\right]D^{j}(\bar{x}_{r}^{N})^{T}=(x^{N}_{r}(t))^{(j)}|_{t=-1}$

Therefore, (3.24) implies (4.25) if the interval between
$\underline{{\boldsymbol{b}}}_{j}$ and ${\boldsymbol{b}}_{j}$ is large enough.
In addition, the spectral coefficients of $\bar{x}^{N}_{r}(D^{T})^{m_{1}}$ is
the same as the spectral coefficients of $(x^{N}_{r}(t))^{(m_{1})}$. From
(3.26) and (3.6), we have

$\displaystyle\sum_{n=0}^{N-r-m_{1}+1}|a^{N}_{n}(m_{1})|\leq{\boldsymbol{d}}$

So, $\\{(\bar{x}^{N},\bar{u}^{N})\\}$ satisfies (4.26). Because we select
$\underline{{\boldsymbol{b}}}$ and $\bar{\boldsymbol{b}}$ large enough so that
the optimal trajectory of the original continuous-time problem is contained in
the interior of the region, we can assume that $(x(t),u(t))$ is also bounded
by $\underline{{\boldsymbol{b}}}$ and $\bar{\boldsymbol{b}}$. Then, (3.23) and
(3.25) imply (4.23) for $N$ large enough. To summarize,
$(\bar{x}^{N},\bar{u}^{N})$ satisfies (4.18)-(4.26). Therefore, it is a
feasible trajectory of (4.18)-(4.23). $\Box$

Proof of Theorem 4: Consider $\\{(\bar{x}^{\ast N},\bar{u}^{\ast
N})\\}_{N=N_{1}}^{\infty}$, a sequence of optimal solutions of Problem ${\rm
B}^{\rm N}$. From Lemma 6,

$\\{||\ddot{x}_{r}^{\ast N}(t)||_{\infty}|N\geq N_{1}\\}$

is bounded. Now, we can apply Lemma 9 to conclude that there exists a
subsequence of $\\{(x^{\ast N}(t),u^{\ast N}(t))\\}_{N=N_{1}}^{\infty}$ and an
optimal solution of Problem B so that the limits in (4.31) converge uniformly.
$\Box$

##  5 Simulation results

The rate of convergence for the optimal cost is illustrated in the following
example

$\displaystyle\min_{u}\int_{0}^{\pi}(1-x_{1}+x_{1}x_{2}+x_{1}u)^{2}dt$ subject
to $\displaystyle\dot{x}_{1}=-x_{1}^{2}x_{2}$
$\displaystyle\dot{x}_{2}=-1+\frac{1}{x_{1}}+x_{2}+\sin t+u$ $\displaystyle
x(0)=\left[\begin{array}[]{cccccccccccccc}1\\\
0\end{array}\right],\;\;x(\pi)=\left[\begin{array}[]{cccccccccccccc}\frac{1}{\pi+1}\\\
2\end{array}\right]$

The analytic solution of this problem is known so that the approximation error
can be computed

$\displaystyle x_{1}(t)=\frac{1}{1-\sin t+t}$ $\displaystyle x_{2}(t)=1-\cos
t$ $\displaystyle u=-(t+1)+\sin t+\cos t$ $\displaystyle\mbox{optimal cost}=0$

The problem is solved using PS optimal control method. The approximated
optimal cost is compared to the true value. The number of nodes, N, ranges
from $4$ to $16$. The error decreases rapidly as shown in Table 1.

N | 4 | 6 | 8 | 10 | 12 | 14 | 16
---|---|---|---|---|---|---|---
Error | $7.5\times 10^{-2}$ | $1.1\times 10^{-3}$ | $2.1\times 10^{-4}$ | $7.1\times 10^{-5}$ | $6.7\times 10^{-6}$ | $1.0\times 10^{-6}$ | $5.8\times 10^{-7}$
Table 1: The error of optimal cost

The rate of convergence is illustrated in the following Figure 3.

Figure 3: Log-scale plot of the error of optimal cost (the solid curve)

Because the analytic solution is $C^{\infty}$, the rate of convergence of the
PS method is faster than any polynomial rate. As a result, it converges
exponentially. Of course, in practical computations the accuracy is limited by
the machine precision of the computers. Therefore, the accuracy cannot be
improved after $N$ is sufficiently large.

##  6 Conclusions

It is proved that the PS optimal control has a high-order rate of convergence.
According to the theorems in Section 3, the approximate cost computed using
the Legendre PS method converges at an order determined by the smoothness of
the original problem. More specifically, the rate is about
$\frac{1}{N^{2m/3-1}}$, where $m$ is defined by the smoothness of the optimal
trajectory. If the cost function can be accurately computed, then the
convergence rate is improved to $\frac{1}{N^{2m-1}}$. If the optimal control
is $C^{\infty}$, then the convergence rate can be made faster than any given
polynomial rate. The results in Section 4 imply that the discretization using
the Legendre PS method is feasible; and there always exists a convergent
subsequence from the approximate discrete optimal solutions, provided some
smoothness assumptions are satisfied.

## References

  * [1] J. T. Betts, Practical Methods for Optimal Control Using Nonlinear Programming, SIAM, Philadelphia, PA, 2001.
  * [2] J. T. Betts, “Survey of Numerical Methods for Trajectory Optimization,” Journal of Guidance, Control, and Dynamics, Vol. 21, No. 2, 1998, pp. 193-207.
  * [3] J. P. Boyd, Chebyshev and Fourier Spectral Methods, second edition, Dover, 2001,
  * [4] C. Canuto, M. Y. Hussaini, A. Quarteroni and T. A. Zang, Spectral Method in Fluid Dynamics. New York: Springer-Verlag, 1988.
  * [5] G. Elnagar, and M. A. Kazemi, Pseudospectral Chebyshev Optimal Control of Constrained Nonlinear Dynamical Systems, Computational Optimization and Applications, 11, 1998, pp. 195-217.
  * [6] Fahroo, F., Ross, I. M., ”Costate Estimation by a Legendre Pseudospectral Method,” Proceedings of the AIAA Guidance, Navigation and Control Conference, 10-12 August 1998, Boston, MA.
  * [7] Q. Gong, W. Kang, and I. M. Ross, A Pseudospectral Method for the Optimal Control of Constrained Feedback Linearizable Systems, IEEE Trans. Automat. Contr., Vol. 51, No. 7, pp. 1115-1129, 2006.
  * [8] Q. Gong, M. Ross, W. Kang, F. Fahroo, Connections Between the Covector Mapping Theorem and Convergence of Pseudospectral Methods for Optimal Control, Computational Optimization and Applications, to appear.
  * [9] W. W. Hager, Runge-Kutta methods in optimal control and the transformed adjoint system, Numerische Mathematik, Vol. 87, pp. 247-282, 2000.
  * [10] A. L. Dontchev and W. W. Hager, The Euler approximation in state constrained optimal control, Mathematics of Computation, Vol. 70, pp. 173-203, 2000.
  * [11] J. Hesthaven, S. Gottlieb, and D. Gottlieb, Spectral Methods for Time-Dependent Problems, Cambridge University Press, 2007.
  * [12] W. Kang, N. Bedrossian, Pseudospectral Optimal Control Theory Makes Debut Flight - Saves NASA $1M in under 3 hrs, SIAM News, September, 2007.
  * [13] W. Kang, Q. Gong, and I. M. Ross, On the Convergence of Nonlinear Optimal Control using Pseudospectral Methods for Feedback Linearizable Systems, International Journal of Robust and Nonlinear Control, Vol. 17, 1251-1277, online publication, 3 January, 2007.
  * [14] W. Kang, I. M. Ross, Q. Gong, Pseudospectral Optimal Control and Its Convergence Theorems, Analysis and Design of Nonlinear Control Systems - In Honor of Alberto Isidori, A. Astolfi and L. Marconi eds., Springer, 2008.
  * [15] S. W. Paris and C. R. Hargraves, OTIS 3.0 Manual, Boeing Space and Defense Group, Seattle, WA, 1996.
  * [16] S. W. Paris, J. P. Riehl, and W. K. Sjauw, “Enhanced Procedures for Direct Trajectory Optimization Using Nonlinear Programming and Implicit Integration,” Proceedings of the AIAA/AAS Astrodynamics Specialist Conference and Exhibit, 21-24 August 2006, Keystone, CO. AIAA Paper No. 2006-6309.
  * [17] E. Polak, Optimization: Algorithms and Consistent Approximations, Springer-Verlag, Heidelberg, 1997.
  * [18] J. P. Riehl, S. W. Paris, and W. K. Sjauw, “Comparision of Implicit Integration Methods for Solving Aerospace Trajectory Optimization Problems,” Proceedings of the AIAA/AAS Astrodynamics Specialist Conference and Exhibit, 21-24 August 2006, Keystone, CO. AIAA Paper No. 2006-6033.
  * [19] Ross, I. M., A Beginner s Guide to DIDO: A MATLAB Application Package for Solving Optimal Control Problems, Elissar Inc., Monterey, CA, October 2007.
  * [20] G. Sansone, A. H. Diamond, and E. Hille, Orthogonal Functions, Robert E. Krieger Publishing Co., Huntington, New York, 1977.
  * [21] A. Wouk, A Course of Applied Functional Analysis, John Wiley & Sons, New York, 1979.