Optimal Control Problem

Authors
Affiliation

Enrico Bertolazzi

University of Trento, Department of Industrial Engineering

Matteo Dalle Vedove

University of Trento, Department of Industrial Engineering

Introduction

The optimal control problems can be seen as a functional minimization, called target, subjected to ordinary differential constraints.

Isoperimetric problem

Let’s consider a function

y (x) \in [a,b]\mapsto\Bbb{R}

where the are known the point y(a) = y(b) = 0. If we consider the function as a rope having fixed length l that, by calculus I, can be computed as

l = \int_a^b \sqrt{1+y'(x)^2}\, \mathrm{d}x

an example of optimal control problem is the one to determine the function f that minimize it’s integral over the it’s domain:

\begin{aligned} \textrm{minimize:} \qquad & \mathcal{F}(y) = - \int_a^b y(x)\, \mathrm{d}x \\ \textrm{subject to:} \qquad & y(a) = 0,\quad y(b) = 0\\ & \int_a^b \sqrt{1+y'(x)^2}\, \mathrm{d}x - l = 0 \end{aligned}

Discretized solution

A way to solve the problem is by discretization, dividing the domain [a,b] in n edges having length h= \frac{b-a}{n} and so x_k = a + k h. By founding the discrete solution y_k of the system we can approximate also the continuous solution as

y_k \approx y(x_k)

Considering the integral approximation

\int_{x_{k-1}}^{x_k} y(x)\, \mathrm{d}x \approx h y_{k-\frac 1 2} = h \frac{y_{k} + y_{k-1}}{2}

and the derivative

y'(x_k) \approx y_{k-\frac 1 2}' = \frac{y_k - y_{k-1}}{2}

With this being stated the initial isoperimetric problem can be rewritten (with some analytical simplification for the calculus) in the discretized form as

\begin{aligned} \textrm{minimize:} \qquad & - \frac h 2\big(y_0+y_n\big) - h \sum_{k=1}^{n-1} y_k \\ \textrm{subject to:} \qquad & y_0 = 0,\quad y_n= 0\\ & l - h \sum_{k=1}^{n} \sqrt{1 + \left( \frac{y_k - y_{k-1}}{h} \right)^2} \end{aligned}

As here stated, this problems become a constrained minimization one that can so be solved using the Lagrange multiplier method on which the unknowns is the vector of the discretized function \bm{y} = (y_1,\dots, y_n) and the 3 Lagrange multipliers \lambda_0, \lambda_n, \lambda associated to the constraints:

\begin{aligned} \mathcal{L}\big(\bm{y},\lambda_0,\lambda_n,\lambda\big) & = - \frac{h}{2}\sum_{k=1}^{n} \big(y_k + y_{k-1}\big) - \lambda_0y_0 - \lambda_ny_n \\ & - \lambda \left( l - h \sum_{k=1}^{n} \sqrt{1 + \left( \frac{y_k - y_{k-1}}{h} \right)^2} \right) \end{aligned}

By calculating the gradient of the lagrangian and setting it to zero, we satisfy the first necessary condition for the solutions that relates to the following non-linear system that can be solved by approximate solutions:

for k=0:

- \dfrac{h}{2} - \lambda_0 + \lambda \dfrac{ \dfrac{y_1-y_0}{h} }{\sqrt{1- \left( \dfrac{y_1-y_0}{h} \right)^2}} = 0

for k = 1,2,\dots,n-1

-h + \dfrac{\lambda \dfrac{y_k-y_{k+1}}{h}}{\sqrt{1+\left( \dfrac{y_k - y_{k-1}}{h} \right)^2}} - \dfrac{\lambda \dfrac{y_{k+1}-y_{k}}{h}}{\sqrt{1+\left( \dfrac{y_{k+1} - y_{k}}{h} \right)^2}} = 0

for k=n:

- \frac{h}{2} - \lambda_n - \lambda \dfrac{\dfrac{y_n - y_0}{h}}{\sqrt{1 + \left( \dfrac{y_1-y_0}{h} \right)^2}} = 0

multipliers conditions:

\begin{cases} y_0 = 0,\qquad y_n = 0 \\ l-h\displaystyle\sum_{k=1}^{n} \sqrt{1 + \left( \dfrac{y_k - y_{k-1}}{h} \right)^2} = 0 \end{cases}

Continuous interpretation

Considering now the limit for h\to 0 the discretized systems seems to become continuous; in particular the 4-th conditions become

y_0 \to y(a) = 0 \quad \textrm{and} \quad y_n \to y(b) = 0

The first condition (associated to k=0) and the third (k = n) instead becomes

- \frac{h}{2} - \lambda_0 + \lambda \dfrac{ \dfrac{y_1-y_0}{h} } {\sqrt{1- \left( \dfrac{y_1-y_0}{h} \right)^2}} \quad\to\quad - \lambda_0 - \lambda \dfrac{y'(a)}{\sqrt{1-y'(a)^2}} = 0

- \frac{h}{2} - \lambda_n - \lambda \dfrac{\dfrac{y_n - y_0}{h}}{\sqrt{1 + \left( \dfrac{y_1-y_0}{h} \right)^2}} \quad\to\quad - \lambda(b) - \lambda \dfrac{y'(b)}{\sqrt{1-y'(b)^2}} = 0

With some analytical manipulation the second lagrange condition can be stated as

- 1 - \lambda \frac{\mathrm{d}}{\mathrm{d}x} \left. \frac{y'(x)}{\sqrt{1+y'(x)^2}} \right|_{x=x_k} = 0

Observation

By pushing the limit h\to 0 the discretized problem becomes continuous and can relate to the general formulation of minimization of a functional \mathcal{F} subjected to an equality constraints (depending on \mathcal G):

\begin{aligned} \textrm{minimize:} \quad & \mathcal{F}(y) = - \int_a^b y(x)\, \mathrm{d}x && \implies\int_a^b L(y,y',t)\,\mathrm{d}t \\ \textrm{subject to:} \quad & y(a) = 0, \quad y(b) = 0\\ & \int_a^b \sqrt{1+y'(x)^2}\, \mathrm{d}x - l = 0 && \implies\int_a^b \mathcal G(y,y',y)\, \mathrm{d}t = 0 \end{aligned}

In the isoperimetric problem we can determine the functions L(y,y',x)=-y whose derivative are

\frac{\partial L}{\partial y} = -1 \quad\textrm{and}\quad \frac{\partial L}{\partial y'} = 0

Using so the Euler-Lagrange approach we can see that it fails because

\frac{\partial L}{\partial y} - \frac{\mathrm{d}}{\mathrm{d}x} \frac{\partial L}{\partial y'} = - 1 \neq 0

Considering instead now the definition of the constraint

\mathcal{G}(y,y',x) = \frac{l}{b-a} - \sqrt{1 + y'^2}

we can build another lagrangian \tilde{ L } that also considers (with a Lagrange multiplier \lambda) the constraint:

\tilde{ L }(y,y',x,\lambda) = L(y,y',x) - \lambda \mathcal{G}(y,y',x)

Applying the Euler-Lagrange method on this expression it’s possible to get the same result achieved with discretization, in fact

\frac{\partial\tilde{L}}{\partial y} - \frac{\mathrm{d}}{\mathrm{d}x} \frac{\partial\tilde{L}}{\partial y'} = - 1 - \lambda \frac{y'}{\sqrt{1 + y'^2}}

General formulation from Euler-Lagrange

Considering the optimal control problem as

\begin{aligned} \textrm{minimize:} \qquad & \int_a^b L (x,x',t)\, \mathrm{d}t \\ \textrm{subject to:}\qquad & x(a) = x_a,\quad x(b) =x_b \\ & \int_a^b \mathcal G(x,x',t)\,\mathrm{d}t = 0 \end{aligned}

a way to reach a solution is by constructing the lagrangian \mathcal{L} defined as

\begin{aligned} \tilde{L}(x,x',t, \lambda) & = L(x,x',t) - \lambda \mathcal G(x,x',t) \\ \mathcal{L}(x,\lambda,\lambda_a,\lambda_b) & = \int_a^b \tilde{ L }(x,x',t,\lambda)\,\mathrm{d}t - \lambda_a x(a) - \lambda_b x(b) \end{aligned}

At this point the solution of the problem is the stationary point (that in reality is a function) x of the lagrangian, and this means setting to zero the Gateaux derivative of \mathcal{L}:

\begin{aligned} \delta\mathcal{L} &= \int_a^b \left( \frac{\partial\tilde{L}}{\partial x}\delta x + \frac{\partial\tilde{L}}{\partial x'} \delta x' + \frac{\partial\tilde{L}}{\partial \lambda} \delta\lambda \right)\, \mathrm{d}t \\ & - \delta\lambda_a \big(x(a)-x_a\big) - \lambda_a \delta x(a) - \delta\lambda_b \big(x(b)-x_b\big) - \lambda_b \delta x(b) \end{aligned}

Performing the integration by part in order to remove the term \delta x' the derivative so becomes

\begin{aligned} \delta \mathcal{L} & = \int_a^b A\,\delta x\, \mathrm{d}x - \delta\lambda_a \big(x(a)-x_a\big) - \delta\lambda_b\big(x(b)-x_b\big) \\ & + B \delta x(a) + C\delta x(b) \end{aligned}

where

A = \frac{\partial\tilde{L}}{\partial x} - \frac{\mathrm{d}}{\mathrm{d}t} \frac{\partial\tilde{L}}{\partial x'} \quad B = -\lambda_a - \frac{\partial\tilde{L}}{\partial y'}\Big|_a \quad C = -\lambda_b + \frac{\partial\tilde{L}}{\partial y'}\Big|_b

At this point the overall solution of the minimization problem becomes the following ordinary differential equation system:

\begin{cases} \left(\dfrac{\partial L}{\partial x} \lambda \dfrac{\partial\mathcal{G}}{\partial x} \right) - \dfrac{\mathrm{d}}{\mathrm{d}x} \left(\dfrac{\partial L}{\partial x'}- \lambda \dfrac{\partial\mathcal{G}}{\partial x'} \right) = 0 \\[1em] x(a) = x_a, \quad x(b) = x_b \\[1em] -\lambda_a - \left( \dfrac{\partial L}{\partial x'}\Big|_a- \lambda \dfrac{\partial\mathcal{G}}{\partial x'} \Big|_a \right) = 0\\[1em] -\lambda_b + \left(\dfrac{\partial L}{\partial x'}\Big|_b - \lambda \dfrac{\partial\mathcal{G}}{\partial x'} \Big|_b \right) = 0 \end{cases}

Ordinary differential equation constraint

Let’s consider the general formulation of the isoperimetric problem as

\begin{aligned} \textrm{minimize:} \qquad & \, \mathcal{F}(x) = \int_a^b L (x,x',t)\, \mathrm{d}t \\ \textrm{subject to:}\qquad & w\big(x(a),x(b) \big)= 0 \\ & \int_a^b \mathcal G(x,x',t)\,\mathrm{d}t = l \end{aligned}

the problem can be solved introducing the function z(t) such that z' = \mathcal{G}(x,x',t): this means that z(a) = 0 and z(b) = l and so the problem if transformed, considering also the introduction of the relation y = x', to the form

\begin{aligned} \textrm{minimize:} \qquad & \mathcal{F}(x,y) = \int_a^b L (x,y,t)\, \mathrm{d}t \\ \textrm{subject to:}\qquad & w\big(x(a),x(b) \big)= 0 \\ & z(a) = 0,\quad z(b) = l \\ & x' = y \\ & z' = \mathcal G(x,y,t) \end{aligned}

With the problem as here stated we refer to x and z as states, while y is the control. We can so see that an optimal control problem can be easily transform into a boundary value problem that can be solved using the technique shown on functional minimization.

General formulation

Considering the problem

\begin{aligned} \textrm{minimize:} \quad & \underbrace{ \overbrace{ \phi\big(\bm{x}(a), \bm{x}(b)\big)}^\textrm{Mayer} + \overbrace{\int_a^b L (\bm{x}, \bm{u},t)\, \mathrm{d}t}^\textrm{ Lagrange}}_\textrm{Bolza} \quad && \textrm{: target} \\ \textrm{subject to:} \quad & \bm{x}' = \bm{f}\big(\bm{x},\bm{u},t\big) && \textrm{: dynamical system} \\ & \bm{b}\big(\bm{x}(a),\bm{x}(b)\big) = 0 && \textrm{: boundary conditions} \\ & \int_a^b \bm{g}(\bm{x},\bm{u}, t)\, \mathrm{d}t = \bm{g}_0 && \textrm{: integral constraints} \end{aligned}

where \bm{x} is the vector of states and \bm{u} the vector of the controls. To solve this kind of problem the first thing to do is to remove the integral constraints transforming them in simple boundary conditions; in order to so is the one to consider g(\bm{x},\bm{u},t) = z' (adding so a new equation in the relation associated to the dynamical system) as a new derivate variable, then by integration it’s possible to see that

z(a) = 0 \qquad z(b) = g_0

With that said the minimization problem becomes

\begin{aligned} \textrm{minimize:} \qquad & \phi\big(\bm{x}(a), \bm{x}(b)\big) + \int_a^b L (\bm{x}, \bm{u},t)\, \mathrm{d}t \\ \textrm{subject to:} \qquad & \bm{x}' = \bm{f}\big(\bm{x},\bm{u},t\big) \\ & \bm{z}' = \bm{g}(\bm{x},\bm{u},t) \\ & \bm{b}\big(\bm{x}(a),\bm{x}(b)\big) = 0 \\ & \bm{z}(a) = 0 \quad \bm{z}(b) = \bm{g}_0 \end{aligned}

Compacting the vector of states \bm{x} and added variable \bm{z} we can define a new variable vector \bm{w}, and similarly all the dynamical system can be described by a multi-variable function F and boundary conditions B and so the problem can be compacted as

\begin{aligned} \textrm{minimize:} \quad & \psi\big(\bm{w}(a), \bm{w}(b)\big) + \int_a^b \mathcal M(\bm{w}, \bm{u},t)\, \mathrm{d}t \\ \textrm{subject to:} \quad & \bm{w}' = \bm{F}\big(\bm{w},\bm{u},t\big) \\ & \bm{B}\big(\bm{w}(a),\bm{w}(b)\big) = 0 \end{aligned}

To solve the problem as here state it’s necessary to build the function of the boundary conditions \bm{\mathcal{B}} (depending by the Lagrange multipliers \bm{\mu}), also known as utility function, and the Hamiltonian \mathcal{H} (depending by the multipliers \bm{\lambda}):

\begin{split} \bm{\mathcal B}\big(\bm{w}(a),\bm{w}(b),\bm{\mu}\big) & = \psi\big(\bm{w}(a),\bm{w}(b)\big) + \bm{\mu}\cdot \bm{B}\big( \bm{w}(a),\bm{w}(b) \big) \\ \mathcal H(\bm{w},\bm{\lambda},\bm{u},t) & = \mathcal M(\bm{w}, \bm{u},t) + \bm{\lambda} \cdot \bm{F}(\bm{w},\bm{u}, t) \end{split}

At this point it’s possible to compute the lagrangian \mathcal{L} on which the variations can be performed:

\begin{aligned} \mathcal{L}\big( \bm{w},\bm{u},\bm{\lambda}, \bm{\mu} ,t \big) & = \bm{\mathcal B}\big(\bm{w}(a),\bm{w}(b),\bm{\mu}\big) + \mathcal H(\bm{w},\bm{\lambda},\bm{u},t) \\ \delta\mathcal{L} & = \int_a^b \Big( A \delta w + B \delta u + C \delta\lambda \Big)\, \mathrm{d}t \\ & + D \, \delta w(a) + E \, \delta w(b) - \bm{B} \big(w(a),w(b)\big) \delta\mu \end{aligned}

where

A = \lambda' + \frac{\partial H}{\partial w} \qquad B = \frac{\partial H}{\partial u} \qquad C = f(w,u,t) - x'

D = \frac{\partial B}{\partial w(a)}\Big|_{a} + \lambda(a) \qquad \qquad E = \frac{\partial B}{\partial w(b)}\Big|_{b} -\lambda(b)

With that said the resulting boundary valued problem becomes

\left\{\begin{aligned} &w' = f(x,u,t) && \textrm{: original ODE} \\ &\lambda' = - \frac{\partial H}{\partial w} && \textrm{: adjoint ODE} \\ &B\big(w(a),w(b)\big) = 0 && \textrm{: original boundary condition} \\ &\left. \begin{aligned} \frac{\partial B}{\partial w(a)}\Big|_{a}+ \lambda(a) = 0\\ \frac{\partial B}{\partial w(b)}\Big|_{b} - \lambda(b) = 0 \end{aligned} \qquad \right\} &&\textrm{: adjoint boundary conditions} \\ & \frac{\partial H}{\partial u} = 0 && \textrm{: control equation} \end{aligned} \right. \tag{1}

Example 1 (optimal control problem) Let’s consider the problem of a mass m sliding on a plane (coordinate x) subjected only to an external applied force F, the optimal control problem is

\begin{aligned} \textrm{minimize:} \qquad & \int_0^1 F^2\, \mathrm{d}t \\ \textrm{subject to:} \qquad & x' = v, \qquad v' = \frac F m \\ & x(0) = 0 \qquad x(1) = 1 \\ & v(0) = 0 \qquad v(1) = 0 \end{aligned}

In this case the description of the dynamical system is known from the physical domain, in fact the velocity is the derivative of the position and the acceleration (derivative of velocity) is equal to the force applied divided by the mass; the boundaries conditions are set by the problem.

At this point to solve the problem we have to compute the two function \mathcal B, \mathcal H defined as

\begin{aligned} \mathcal B & = \mu_1 \big(x(a) - 0\big) + \mu_2\big(x(b)-1\big) + \mu_3\big(v(a)-0\big) + \mu_4 \big(v(b) - 0\big) \\ \mathcal H & = F^2 + \lambda_1 v + \lambda_2 \frac F m \end{aligned}

Following the results of equation Equation 1 the boundary valued problem that has to be solved to find the minimum point is

\begin{cases} x' = v \\ v' = F/m \\ \lambda_1' = 0 \\ \lambda_2' = \lambda_1 \\ x(0) = 0 \qquad x(1) = 1 \\ v(0) = 0 \qquad v(1) = 0 \\ \mu_1 + \lambda_1(0) = 0 \\ \mu_3 + \lambda_2(0) = 0 \\ \mu_2 - \lambda_1(1) = 0 \\ \mu_4 - \lambda_2(1) = 0 \\ 2 F + \frac{\lambda_2}{m} = 0 \end{cases}

In this case the adjoin boundary conditions are trivially solved (in fact will exists \mu_i that will satisfy the equation); considering the last equation in the expression of the velocity the system becomes

\begin{cases} x' = v \\ v' = -\frac{\lambda_2}{2m^2} \\ \lambda_1' = 0 \\ \lambda_2' = \lambda_1 \\ x(0) = 0 \qquad x(1) = 1 \\ v(0) = 0 \qquad v(1) = 0 \\ \end{cases}

Considering that \lambda_1 is a constant c_1 (in order to have null derivative), then it means that \lambda_2 is in the form c_1 t + c_2 and so the expression of the acceleration and velocity becomes

\begin{aligned} v' & = -\frac{c_1 t + c_2}{2m^2} \\ \xrightarrow{\int} \quad v & = - \frac{c_1}{4m^2} t^2 - \frac{c_2}{2m^2}t + c_3 \\ \xrightarrow{\int} \qquad x & = - \frac{c_1}{12m^2}t^3 - \frac{c_2}{4m^2}t^2 + c_3 t + c_4 \end{aligned}

Considering the boundary conditions than the integration constant are c_3 = c_4 = 0 and the solution of the linear system given by -c_1 - 2c_2 = 0 and -c_1 -3c_2 = 12m^2 that determines the last full solution

\begin{aligned} x & = - 6t^2 + 6 t \\ v & = - 2t^3 + 3 t^2 \\ \lambda_1 & = 24m^2 \\ \lambda_2 & = 24m^2 t - 12m^2 \end{aligned}

and so the control history to determine this minimal solution is

F = - \frac{\lambda_2}{2m} = 6\big(1-2t\big) m

Pontryagin maximum (minimum) principle

The Pontryagin maximum (minimum) principle has been developed in order to solve problem where the controls must stay bounded to a certain range due to the physical implementation of the system (if, as example, a motor is powered with 20W of power, than it cannot move object that require more than that power).

Mathematically this kind of problem is described by

\begin{aligned} \textrm{minimize:} \qquad & \phi\big(x(a),x(B)\big) + \int_a^b L (x,u,t) \, \mathrm{d}t \\ \textrm{subject to:} \qquad & x' = f(x,u,t) \\ & b\big(x(a), x(b)\big) = 0 \\ & u(t) \in \mathcal U \end{aligned}

where \mathcal U is the domain of the controls that, in order for the principle to work, must be convex and compact.

To solve this problem we compute the Hamiltonian \mathcal H and the utility function \mathcal{B} as in the previous cases:

\begin{aligned} \mathcal{H}(x,u, \lambda,t) & = L(x,u,t) + \lambda f(x,u,t) \qquad \\ \mathcal{B}(x_a,x_b,\mu) & = \phi(x_a,x_b) + \mu b(x_a,x_b) \end{aligned}

With that stated the boundary valued problem becomes similar to the formulation in Equation 1 with the addition of the Pontryagin minimum principle:

\left\{ \begin{aligned} & x' = f(x,u,t) = \frac{\partial\mathcal H}{\partial\lambda} && \textrm{: original ODE} \\ &\lambda' = - \frac{\partial\mathcal H}{\partial x} && \textrm{: adjoint ODE - co-equations} \\ &b\big(w(a),w(b)\big) = 0 && \textrm{: original BC} \\ &\left. \begin{aligned} \frac{\partial\mathcal{B}}{x_a} + \lambda(a) = \bm{0} \\ \frac{\partial\mathcal{B}}{x_b} - \lambda(b) = \bm{0} \end{aligned} \qquad \right\} &&\textrm{: adjoint BC} \\ & u(t) = \underset{\overline u \in \mathcal U}{\textrm{argmin}} \big\{ \mathcal{H}\big(x, \overline u,\lambda, t\big) \big\} && \textrm{: Pontryagin min principle} \\ & \frac{\partial\mathcal{H}}{\partial u} = 0 && \textrm{: control equation} \end{aligned} \right.

Examples with Pontryagin maximum (minimum) principle

Example 2 (maximum travel optimal control problem) Let’s consider the system in example Example 1 where in this case the goal is to maximize the travel of the mass m having a force F \leq |1|; this means solving the following optimal control problem

\begin{aligned} \textrm{minimize:} \qquad & -x(1) \\ \textrm{subject to:} \qquad & x' = v, \qquad v' = \frac F m \\ & x(0) = 0 \\ & v(0) = 0 \qquad v(1) = 0 \\ & F \in [-1,1] \end{aligned}

In order to solve this problem it’s necessary to determine the Hamiltonian

\mathcal{H}(x,v,\lambda_1,\lambda_2, F,t) = \lambda_1 v + \lambda_2 \frac F m

and the utility function

\mathcal{B}(x_a,v_a,x_b,v_b,\mu_1,\mu_2,\mu_3) = -x_b + \mu_1 x_a + \mu_2 v_a + \mu_3 v_b

The adjoint equation so becomes

\lambda_1' = - \frac{\partial\mathcal{H}}{\partial x} = 0

and

\lambda_2' = - \frac{\partial\mathcal{H}}{\partial v} = - \lambda_1

for the the adjoint boundary conditions only -1 - \lambda_1(1)=0 is reported (the other equation trivilly solves the cofficients \mu_i that gives no information for the final solution). The boundary valued problem so become

\left\{ \begin{aligned} & x' = v \\ & v' = F/m \\ & \lambda_1' = 0 \\ & \lambda_2' = - \lambda_1 \\ & -1-\lambda_1(1)= 0\\ & x(0) = v(0) = v(1) = 0 \\ & F(t) = \underset{\overline F \in [-1,1]}{\textrm{argmin}} \left\{ \lambda_1(t) v + \lambda_2(t) \frac {\overline F}m \right\} \end{aligned} \right.

In this case in order to compute the argument minimum associated to the Pontryagin principle we have to consider only the term related to \overline F, and in particular this means

\begin{aligned} F(t) & = \underset{\overline F \in [-1,1]}{\textrm{argmin}} \left\{ \lambda_2(t) \frac {\overline F}m \right\} \\& = \begin{cases} + 1 & \textrm{if } \lambda_2(t) < 0 \\ - 1 & \textrm{if } \lambda_2(t) > 0 \\ [-1,1] \qquad & \textrm{if } \lambda_2(t) = 0 \\ \end{cases} \\ & = - \textrm{sign} \lambda_2(t) \end{aligned}

and so the simplified system that determines the solution becomes

\left\{ \begin{aligned} & x' = v \\ & v' = -\frac{\textrm{sign}\lambda_2}{m} \\ & \lambda_1' = 0 \\ & \lambda_2' = - \lambda_1 \\ & \lambda_1(1)= -1\\ & x(0) = v(0) = v(1) = 0 \end{aligned} \right.

By integration it can be seen that \lambda_1 can be regarded as a constant c_1 that, for the adjoint boundary condition, is equal to -1 and so \lambda_2 is in the form

\lambda_2 = t + c_2

From now on the solution becomes analytically complex, however it can be found that c_2 = -\frac{1}{2} determining the solution

F(t) = \begin{cases} 1 \qquad & t \leq \frac 1 2 \\ -1 & t > \frac 1 2 \end{cases}


Example 3 (optimal control problem) Given the optimal control problem

\begin{aligned} \textrm{minimize:} \qquad & \mathcal{F} (x,u) = \int_0^2 2x- u^2 \, \mathrm{d}t \\ \textrm{subject to:} \qquad & x' = 1-2u \\ & x(0) = 1 \qquad x(2) = 0 \end{aligned}

the solution can be obtained by computing the hamiltonian \mathcal{H} and the utility function \mathcal{B} related to the problem:

\begin{aligned} \mathcal{H}(x,u,\lambda) & = 2x - u^2 - \lambda \big(1-2u\big) \\ \mathcal{B}(x_0,x_2, \mu_1, \mu_2) & = - \mu_1 \big(x_0-1) - \mu_2 x_2 \end{aligned}

With that stated, the lonely resulting co-equation is \lambda' = \frac{\partial\mathcal{H}}{\partial x} = 2, determining so that the solution of \lambda(t) is in the form 2t + c_1 (where c_1 is a constant depending from the boundary conditions). The adjoint boundary conditions aren’t useful to determine the solution because they trivially solves the parameters \mu_1,\mu_2. We still have the control equation that states

\frac{\partial\mathcal{H}}{\partial u} = -2u + 2\lambda = 0 \qquad \implies\quad u(t) = \lambda(t)

Using the original ordinary differential equation we have that x' = 1 - 4t - 2c_1 that by integration becomes

x(t) = -2t^2 + \big(1-2c_1\big) t + c_2

From the boundary condition x(0) = 1 we obtain c_2 = 1, while considering x(2) = 0 we have c_1 = -5/4 and so the final control law is

u(t) = 2t - \frac{5}{4}


Example 4 (optimal control problem) To solve the optimal control problem

\begin{aligned} \textrm{minimize:} \qquad & \int_0^1 t + x - y^2 \, \mathrm{d}t \\ \textrm{subject to:} \qquad & x' = x+y \qquad y' = x + yu \\ & \int_0^1\sin y\, \mathrm{d}t = 3\\ & x(0) = 1 \qquad y(1) = 0 \\ & u \in [-2,3] \end{aligned}

the first thing to do is removing the integral constraint by transforming him in a differential equation with 2 boundary condition as follows:

\begin{aligned} \textrm{minimize:} \qquad & \int_0^1 t + x - y^2 \, \mathrm{d}t \\ \textrm{subject to:} \qquad & x' = x+y \\ & y' = x + yu \\ & z' = \sin y \\ & x(0) = 1, \qquad & z(0) = 0 \\ & y(1) = 0, \qquad & z(1) = 3 \\ & u \in [-2,3] \end{aligned}

The associated hamiltonian \mathcal{H} and utility function \mathcal{B} of the problem are

\begin{aligned} \mathcal{H}(x,y,z,\lambda_1, \lambda_2, \lambda_3) & = t + x - y^2 \\ & - \lambda_1 (x+y) - \lambda_2 (x+yu) - \lambda_3 (\sin y) \\ \mathcal{B}(x_0, y_1, z_0, z_1, \mu_{1\dots 4}) & = -\mu_1 \big(x_0 - 1\big) - \mu_2 y_1 \\ &- \mu_3 z_0 - \mu_4 \big(z_1 -3\big) \end{aligned}

With that stated the co-equations are computed as

\begin{aligned} \lambda_1' & = - \frac{\partial\mathcal{H}}{\partial x} = -1 + \lambda_1 + \lambda_2 \\ \lambda_2' & = - \frac{\partial\mathcal{H}}{\partial y} = 2y + \lambda_1 + \lambda_2 u + \lambda_3 \cos y \\ \lambda_3' & = - \frac{\partial\mathcal{H}}{\partial z} = 0 \end{aligned}

Evaluating the adjoint boundary conditions for the point that present no explicit constraints (so x_1, y_0) we retrieve \lambda_1(1) = 0 and \lambda_2(0) = 0. With the definition of the control equation \frac{\partial\mathcal{H}}{\partial u} = - \lambda_2 y the resulting boundary value problem of the optimal control one is

\begin{cases} x' = x+y \\ y' = x + yu \\ z' = \sin y \\ \lambda_1' = \lambda_1 + \lambda_2 - 1 \\ \lambda_2' = 2y + \lambda_1 + \lambda_2 y + \lambda_3 \cos y \\ \lambda_3' = 0 \\ x(0) = 1 \qquad \lambda_1(1) = 0 \\ y(1) = 0 \qquad \lambda_2(0) = 0 \\ z(0) = 0 \qquad z(1) = 3 \\ -\lambda_2 = 0 \\ u = \underset{\overline u \in [-2,3]}{\textrm{argmin}} \left\{ -\lambda_2 y \overline u \right\} \end{cases}

where the solution of the Pontryagin minimum principle is

u(t) = \begin{cases} -2 \qquad & \lambda_2 y < 0 \\ 3 & \lambda_2 y \geq 0 \end{cases}


Example 5 (optimal control problem) Given the optimal control problem

\begin{aligned} \textrm{minimize:} \qquad & \mathcal{F} (x,y,u) = x(0) + \int_0^1 x^2 + y^2 + u \, \mathrm{d}t \\ \textrm{subject to:} \qquad & x' = y - u, \qquad y' = x+u \\ &x(1) = 2 \\ & u(t) \in [-1,2] \end{aligned}

the solution can be obtained by firstly computing the hamiltonian \mathcal{H} and the utility function \mathcal{B} of the problems resulting in

\begin{aligned} \mathcal{H}\big(x,y,u,\lambda_1,\lambda_2\big) & = x^2 + y^2 + u + \lambda_1 \big(y-u\big) + \lambda_2 \big(x + u\big) \\ & = x^2 + y^2 + u \big(1 - \lambda_1 + \lambda_2\big) + \lambda_1 y + \lambda_2 x \\ \mathcal{B}\big(x_0,x_1, \mu\big) & = x_0 + \mu \big( x_1 - 2\big) \end{aligned}

We can so explicit: the adjoint ordinary differential equations

\begin{aligned} \lambda_1' &= -\frac{\partial\mathcal{H}}{\partial x} = - 2x - \lambda_2 \\ \lambda_2' &= - \frac{\partial\mathcal{H}}{\partial y} = -2y - \lambda_1 \end{aligned}

the adjoint boundary conditions

\begin{aligned} \lambda_1(0) & = - \frac{\partial\mathcal{B}}{\partial x_0} = - 1 \\ \lambda_1(1) &= \frac{\partial\mathcal{B}}{\partial x_1} = \mu \\ \lambda_2(0) & = -\frac{\partial\mathcal{B}}{\partial y_0} = 0 \\ \lambda_2(1) & = \frac{\partial\mathcal{B}}{\partial y_1} = 0 \end{aligned}

and the control equation

\frac{\partial\mathcal{H}}{\partial u} = 1 - \lambda_1 + \lambda_2 = 0

With the addition of the original differential equation the boundary value problem becomes

\begin{cases} x' = y - u \\ y' = x+u \\ \lambda_1' = -2x - \lambda_2 \\ \lambda_2' = -2y - \lambda_1 \\ \lambda_1(0) = -1 \\ \cancel{\lambda_1(1) = \mu} \qquad \qquad \textrm{: trivially satisfied} \\ \lambda_2(0) = \lambda_2(1) = 0 \\ 1 + \lambda_1 - \lambda_2 = 0 \\ u = \underset{\overline u \in [-1,2]}{\textrm{argmin}} \left\{ \overline u \big(1 + \lambda_1 - \lambda_2\big) \right\} \end{cases}

In this case the solution of the Pontryagin minimum principle is reduced to the form

u(t) = \begin{cases} -1 \quad & 1 + \lambda_1 - \lambda_2 > 0 \\ 2 & 1 + \lambda_1 - \lambda_2 < 0 \end{cases}

Problems with linear controls and singular arc

In this section, we examine problems in which the Hamiltonian is a linear function of the control variable u. Although linear problems are generally expected to be easier to handle, this is not always the case. The complexity arises from the fact that the necessary condition for minimizing the Hamiltonian with respect to the control along an optimal trajectory does not yield a well-defined expression for synthesizing the optimal control.

Another challenge associated with this class of problems is the occurrence of discontinuities in the control, which can be difficult to identify. Additionally, there may be an infinite number of switching points that can accumulate at a specific time instant. Given that general results and theorems regarding existence and uniqueness are limited in this context, it becomes necessary to introduce additional relationships by manipulating other necessary conditions.

Consequently, problems with a Hamiltonian that is linear in terms of the control u represent an almost independent area of research within optimal control theory. Most theoretical results in this domain are derived using a geometric approach, involving fiber bundles and symmetry groups rather than traditional variational techniques.

To maintain clarity and simplicity in notation, we will focus on problems with a single control variable and fixed endpoints. The formulation of these problems is as follows:

Definition 1 Minimize the cost functional J subject to the dynamic \bm{x}' with the control u constrained in |u|\leq 1 and the initial state \bm{x}(0)=\bm{x}_0, where

\begin{aligned} J(u)&= \int_0^T f_0(\bm{x})+b_0(\bm{x})u \,\mathrm{d}t\\ \bm{x}'(t) &= \bm{f}(\bm{x})+\bm{b}(\bm{x})u\qquad \mathrm{i.e.} \end{aligned} \tag{2}

f_0(\bm{x}) and b_0(\bm{x}) are scalar functions of the state \bm{x}, u is the scalar control.

To further simplify the computations, we consider the scalar equation of the functional as an extra state variable posing

x_0'= f_0(\bm{x})+b_0(\bm{x})u, \qquad x_0(0)=0,

(we have already seen that this is not a loss of generality), but the key point here is that the control appears only linearly in both dynamics and cost functional. We can write the Hamiltonian for this system, namely

\mathcal{H}(\bm{x},\bm{\lambda},t) = \sum_{i=0}^n f_i(\bm{x})\lambda_i+u\sum_{i=0}^n b_i(\bm{x})\lambda_i.

It is convenient to rename the two addends as

\begin{aligned} A(\bm{x},\bm{\lambda},t)&=\sum_{i=0}^n f_i(\bm{x})\lambda_i,\\ B(\bm{x},\bm{\lambda},t)&=\sum_{i=0}^n b_i(\bm{x})\lambda_i, \end{aligned}

so that the Hamiltonian becomes

\mathcal{H}(\bm{x},\bm{\lambda},t) = A(\bm{x},\bm{\lambda},t)+B(\bm{x},\bm{\lambda},t)u. \tag{3}

The adjoint variables of the costate are given by the usual relation

\lambda_i' = -\dfrac{\partial\mathcal{H}}{\partial x_i}= -\sum_{j=0}^n\lambda_j\dfrac{\partial f_j(\bm{x})}{\partial x_i}-u\sum_{j=0}^n\lambda_j\dfrac{\partial b_j(\bm{x})}{\partial x_i}.

We state now the necessary conditions for problem Equation 2.

Theorem 1 If u^\star(t) is an optimal control and if \bm{x}^\star(t) is the corresponding optimal trajectory, then there are a costate \lambda_i^\star(t) such that

  • for i=1,\ldots,n \begin{aligned} x_{i}'^{\star} &= f_i(\bm{x}^{\star})+b_i(\bm{x}^{\star})u^{\star},\\ \lambda_i'^{\star} &= -\sum_{j=0}^n\lambda_j^\star\dfrac{\partial f_j(\bm{x}^\star)}{\partial x_i^\star}- u^\star\sum_{j=0}^n\lambda_j^\star(t)\dfrac{\partial b_j(\bm{x}^\star)}{\partial x_{i}^\star}; \end{aligned}

  • for all t\in[0,T] and all u(t) s.t. |u(t)|\leq 1, holds A(\bm{x}^\star,\lambda^\star,t)+u^\star(t)B(\bm{x}^\star, \lambda^\star,t)\leq A(\bm{x}^\star,\lambda^\star,t)+u(t)B(\bm{x}^\star,\lambda^\star,t);

  • If T is free, then for all t\in[0,T]

    \mathcal{H}^\star(\bm{x}^\star, \bm{\lambda}^\star,t) := A(\bm{x}^\star, \bm{\lambda}^\star,t)+u^\star(t)B(\bm{x}^\star, \bm{\lambda}^\star,t)=0 , \tag{4}

    while, if T is fixed, then

    \mathcal{H}^\star(\bm{x}^\star, \bm{\lambda}^\star,t) := A(\bm{x}^\star, \bm{\lambda}^\star,t) +u^\star(t)B(\bm{x}^\star, \bm{\lambda}^\star,t)=c , \tag{5} for a real constant c.

It follows easily from equation Equation 3 with the PMP, that the optimal control is given by

u^\star(t)=-\mathrm{sign}\{B(\bm{x}^\star, \bm{\lambda}^\star,t)\}=-\mathrm{sign}\left\{\sum_{i=0}^n b_i(\bm{x}^\star(t))\lambda^\star_i(t)\right\},

as long as the function B(\bm{x}^\star,\bm{\lambda}^\star,t) is not zero. If B becomes zero, then the sign function is not defined. If B is zero only in a point we speak of bang-bang controls, if B is identically zero for an entire interval t\in(t_1,t_2], the controls are called singular. In this case, the associated trajectory \bm{x}(t) is called a singular arc.

Singular controls

If problem Equation 2 results singular, controls, costate and trajectory have the following property: there is at least one half-open interval (t_1,t_2]\subset[0,T] such that

\sum_{i=0}^n \lambda^\star_i(t)b_i(\bm{x}^\star(t))=0, \qquad t\in(t_1,t_2]. \tag{6}

As before, the existence of an extremal singular control, does not imply that the optimal control is singular: we need additional information as uniqueness to conclude its optimality. The function Equation 6 is often called switching function.

Test for singular controls

Suppose we are in the free terminal time case, so that Equation 4 holds, we test if it is possible to have a singular control as follows. First we need to assume that the switching function is zero in an unknown interval (t_1,t_2], then, because of the free end time hypothesis, for t\in(t_1,t_2],

\begin{aligned} \mathcal{H}(\bm{x},\bm{\lambda},t) = \sum_{i=0}^n f_i(\bm{x})\lambda_i+u\sum_{i=0}^n b_i(\bm{x})\lambda_i &= 0\\ & \Downarrow\\ \sum_{i=0}^n f_i(\bm{x})\lambda_i&=0. \end{aligned}

Moreover, for each k\in\mathbb{N}, we have that

\dfrac{\mathrm{d}^k}{\mathrm{d}t^k} \sum_{i=0}^n \lambda_i(t)b_i(\bm{x})=0, \quad \left\{\begin{aligned} k & =1,2,3,\ldots, \\ t & \in(t_1,t_2], \end{aligned}\right. \tag{7}

and similarly,

\dfrac{\mathrm{d}^k}{\mathrm{d}t^k}\sum_{i=0}^n \lambda_i f_i(\bm{x})=0, \quad \left\{\begin{aligned} k & =1,2,3,\ldots, \\ t & \in(t_1,t_2], \end{aligned}\right. \tag{8}

However, the canonical equations (we omit explicit time dependence) are

\begin{aligned} x_i' &= f_i+ub_i\\ \lambda_i' &= -\sum_{j=0}^n \lambda_j(t)\dfrac{\partial f_j}{\partial x_i} -u \sum_{j=0}^n \lambda_j \dfrac{\partial b_j}{\partial x_i}, \quad \left\{\begin{aligned} i,j & =1,2,\ldots,n \\ k & \in\Bbb{N}, \end{aligned}\right. \end{aligned}

Let k=1 in Equation 7, then by the chain rule we have

\dfrac{\mathrm{d}}{\mathrm{d}t}\sum_{i=0}^n \lambda_i b_i = \sum_{i=0}^n\left(\lambda_i'b_i+\lambda_i \sum_{j=0}^n \dfrac{\partial b_j}{\partial x_i}x'_j\right)=0.

Substituting the canonical equations in the previous one, yields

\sum_{i=0}^n\sum_{j=0}^n\left( f_j\lambda_i\dfrac{\partial b_i}{\partial x_j}-b_i\lambda_j\dfrac{\partial f_j}{\partial x_i}\right) +u \sum_{i=0}^n\sum_{j=0}^n\left( \lambda_i b_j\dfrac{\partial b_i}{\partial x_j}- \lambda_j b_i\dfrac{\partial b_j}{\partial x_i}\right)=0.

We notice that the coefficient of u is zero, and we conclude that

\sum_{i=0}^n\sum_{j=0}^n\lambda_i \left( b_j\dfrac{\partial f_i}{\partial x_j}-f_j\dfrac{\partial b_i}{\partial x_j}\right)=0. \tag{9}

Applying the same argument to equation Equation 8, gives

\begin{aligned} 0 &= \dfrac{\mathrm{d}}{\mathrm{d}t}\sum_{i=0}^n \lambda_if_i(\bm{x}) \\ &= \sum_{i=0}^n\left(\lambda_i'f_i+\lambda_i\sum_{j=0}^n\dfrac{\partial f_i}{\partial x_j}x_j'\right)\\ &= \sum_{i=0}^n\sum_{j=0}^n\left( \lambda_if_i\dfrac{\partial f_i}{\partial x_j}-\lambda_jf_i\dfrac{\partial f_j}{\partial x_i}\right) +u \sum_{i=0}^n\sum_{j=0}^n\left( \lambda_i b_j\dfrac{\partial f_i}{\partial x_j}- \lambda_j f_i\dfrac{\partial b_j}{\partial x_i}\right)\\ &= u \sum_{i=0}^n\sum_{j=0}^n\lambda_i\left( b_j\dfrac{\partial f_i}{\partial x_j}-f_j\dfrac{\partial b_i}{\partial x_j}\right)=0. \end{aligned}

This last line implies that either u=0 or

\sum_{i=0}^n\sum_{j=0}^n\lambda_i\left( b_j\dfrac{\partial f_i}{\partial x_j}-f_j\dfrac{\partial b_i}{\partial x_j}\right)=0.

Noticing that the previous line is equal to equation Equation 9, we can not conclude that u=0, and this is true for all k\in\mathbb{N}.\ In conclusion, a necessary but not sufficient test for an extremal singular control is that in the interval (t_1,t_2] the following relations are satisfied: \left.\begin{aligned} \displaystyle\sum_{i=0}^n b_i\lambda_i=0\\ \displaystyle\sum_{i=0}^n f_i\lambda_i=0\\ \displaystyle\sum_{i=0}^n\sum_{j=0}^n\lambda_i\left( b_j\dfrac{\partial f_i}{\partial x_j}-f_j\dfrac{\partial b_i}{\partial x_j}\right)=0 \end{aligned}\right\}

We show now how to compute the explicit expression for the singular control, we return back to the equation of the switching function Equation 7 and let k=2,3,\ldots. It is shown in Athans and Falb (1966, pag.487) that those derivatives for k=2,3,\ldots require extensive manipulations but have all the following same structure:

\sum_{i=0}^n\lambda_i\psi_{ki}(\bm{x})+u\sum_{i=0}^n\lambda_i\phi_{ki}(\bm{x})=0.

If \sum \lambda_i\phi_{ki}(\bm{x})=0 we can perform successive derivatives of Equation 7, in general there will be a finite value of k=m for which \sum \lambda_i\phi_{mi}(\bm{x})\neq 0 and then the singular control u can be solved:

u= -\dfrac{\displaystyle\sum_{i=0}^n\lambda_i\psi_{mi}(\bm{x})}{\displaystyle\sum_{i=0}^n\lambda_i\phi_{mi}(\bm{x})}.

Up to now we have dealt with free final time, if the final time is fixed the we have to take in account the presence of the constant c in the Hamiltonian, but Equation 7 and Equation 8 still hold and hence also the conclusions.\ In practical problems, we must check if

\sum \lambda_i\phi_{ki}(\bm{x})=0

for all k=1,\ldots,m-1 and

\sum \lambda_i\phi_{mi}(\bm{x})\neq 0

If any of the relations are violated then this represent a violation of the necessary conditions and so singular controls can not occur, if these conditions are all met, then there may be a singular extremal control.

References

Athans, Michael, and Peter L. Falb. 1966. Optimal Control : An Introduction to the Theory and Its Applications. Lincoln Laboratory Publications. New York, Saint Louis, San Francisco: McGraw-Hill. http://opac.inria.fr/record=b1080535.
Betts, John T. 2010. Practical Methods for Optimal Control Using Nonlinear Programming. 3rd ed. Society for Industrial; Applied Mathematics.
Bryson, Arthur E., and Yu-Chi Ho. 1975. Applied Optimal Control: Optimization, Estimation, and Control. Wiley.
Kirk, Donald E. 2004. Optimal Control Theory: An Introduction. Dover Publications.
Liberzon, Daniel. 2012. Calculus of Variations and Optimal Control Theory: A Concise Introduction. Princeton University Press.