# Convergence of machine learning methods for feedback control laws:
averaged feedback learning scheme and data driven method

Karl Kunisch  Donato Vásquez-Varas

###### Abstract

This work addresses the synthesis of optimal feedback control laws via machine
learning. In particular, the _Averaged Feedback Learning Scheme_ (AFLS) and a
data driven method are considered. Hypotheses for each method ensuring the
convergence of the evaluation of the objective function of the underlying
control problem at the obtained feedback-laws towards the optimal value
function are provided. These hypotheses are connected to the regularity of the
value function and the stability of the dynamics. In the case of AFLS these
hypotheses only require Hölder continuity of the value function, whereas for
the data driven method the value function must be at least $C^{2}$. It is
demonstrated that these methods are connected via their optimality conditions.
Additionally, numerical experiments are provided by applying both methods to a
family control problems, parameterized by a positive real number which
controls the regularity of the value function. For small parameters the value
function is smooth and in contrast for large parameters it is non-
differentiable, but semi-concave. The results of the experiments indicate that
both methods have a similar performance for the case that the value function
is smooth. On the other hand, if the value function is not differentiable,
AFLS has a better performance which is consistent with the obtained
convergence results.

###### keywords:

optimal feedback control, Hamilton-Jacobi-Bellman equation, learning approach,
data driven.

††journal: Arxiv

[RICAM]organization=Radon Institute for Computational and Applied Mathematics,
Austrian Academy of Sciences, addressline=Altenbergerstraße 69, city=Linz,
postcode=4040, state=Upper Austria, country=Austria
[Graz]organization=Institute of Mathematics and Scientific Computing,
University of Graz, addressline=Heinrichstraße 36, city=Graz, postcode=8010,
state=Styria, country=Austria

AMS classification: 49L12, 49J15, 49N35, 68Q32, 93B52,

##  1 Introduction

The problem of constructing feedback-laws via machine learning methods has
attracted notable attention. This is due to the fact that classical method for
synthesizing feedback laws suffer from the _curse of dimensionality_. Among
the well-known methods for solving the HJB equation we mention finite
difference schemes [1], semi-Lagrangian schemes [2], and policy iteration [3,
4, 5, 6]. The curse of dimensionality can be elevated by machine learning
methods, since they are capable of efficiently approximating high dimensional
functions. For instance we can mention the following contributions:
representation formulas [7, 8, 9, 10], approximating the HJB equation by
neural networks [11, 12, 13, 14, 15, 16, 17, 18, 19], data driven approaches
[20, 21, 22, 23, 24, 25], max-plus methods [26, 27, 28], polynomial
approximation [29, 30], tensor decomposition methods [31, 32, 33, 34, 35, 36],
POD methods [37, 38], tree structure algorithms [39], and sparse grids
techniques[40, 41, 42, 43], see also the proceedings volume [44]. Frequently,
there is no proof concerning the question whether the feedback law constructed
via machine learning provides an optimal control, unless the value function is
sufficiently smooth as is the case in [25], where the value functions is
supposed to be an element of a reproducing kernel Hilbert space.

In the present work we analyze and study the convergence of two methods and
investigate the interconnections between them. The first one consists in
finding a feedback law minimizing the average of the cost of the control
problem with respect to a set of initial conditions. We call this method
_Averaged Feedback Learning Scheme_ (AFLS). It has been introduced and used in
earlier work, see e.g. [16, 45, 45, 46, 14, 47, 17], but has not been
attributed a name before. The second one approximates the feedback law by
using a regression or interpolation algorithm. In each case the feedback
control is parameterized by means of the verification theorem.

The convergence of the methods is investigated in terms of the distance
between the evaluation of the objective function obtained by the feedback law
of the respective method and the value function of the underlying optimal
control problem. This is explained in more detail this in Section 3. We are
able to prove the convergence of the first method under the hypothesis of the
existence of a stable sequence of consistent feedback laws , see[46]. This in
turn depends on the regularity of the value function. For the second method we
are able to prove its convergence if the value function is $C^{2}$.

In order to illustrate the dependence of the convergence in dependence of the
regularity of the value function, we also study numerically the convergence of
the different methods for a family of control problem indexed by a parameter
$\gamma>0$. This problem was first given in [46] where it was proved that its
value function is Lipschitz continuous for all $\gamma>0$ and is non-
differentiable for $\gamma$ large enough. Additionally, in the present article
we prove that for $\gamma$ small enough the value function is
$C^{2}(\mathbb{R}^{d})$. Thus applying the methods for this problem for
different values of $\gamma$ allows us to observe the influence of the
regularity on the behavior of the two proposed methods.

##  2 Control problem

We consider the task of finding an approximate optimal feedback-law for the
following control problem

$\min_{u\in
L^{2}((0,\infty);\mathbb{R}^{m})}\int_{0}^{\infty}\ell(y(t;y_{0},u))+\frac{\beta}{2}|u|^{2}dt,$
(2.1)

where $y(\cdot;y_{0},u)$ is the solution of the equation

$y^{\prime}=f(y)+Bu,\ y(0)=y_{0}.$ (2.2)

Throughout we assue that $B\in\mathbb{R}^{d\times m}$, that $f$ and $\ell$ are
Lipschitz continuous on bounded subsets of $\mathbb{R}^{d}$, and
$\ell\geqslant 0$.

A classical way to find a feedback-law is based on dynamic programming. For
this purpose the value function associated to (2.1) - (2.2) is utilized. It is
defined by

$V(y_{0})=\min_{u\in
L^{2}((0,\infty;\mathbb{R}^{d}))}\int_{0}^{\infty}\left(\ell(y(t;y_{0},u))+\frac{\beta}{2}|u(t)|^{2}\right)dt.$
(2.3)

In the case that $V$ is $C^{1,1}(\mathbb{R}^{d})$ it is known that the
feedback law given by

$u(y)=-\frac{1}{\beta}B^{\top}\nabla V(y)\mbox{ for }y\in\mathbb{R}^{d}$

is optimal, i.e., for $y_{0}\in\mathbb{R}^{d}$, denoting by $y^{*}\in
H^{1}_{loc}((0,\infty);\mathbb{R}^{d})$ the unique solution of

$y^{\prime}(t)=f(y(t))+Bu(y(t)),\ t>0,\ y(0)=y_{0}$

and setting $u^{*}(t)=u(y^{*}(t))$ for $t>0$, we have that $u^{*}$ is an
optimal solution of (2.1). Further, in this case $V$ is a classical solution
of the following Hamilton-Jacobi equation

$\max_{u\in\mathbb{R}^{m}}\left\\{-\nabla
V(y)(f(y)+Bu)-\ell(y)-\frac{\beta}{2}|u|^{2}\right\\}=0\mbox{ for all
}y\in\mathbb{R}^{d}.$ (2.4)

Hence, for finding an optimal feedback law it is enough to solve (2.4). But
there are two main obstacles for using this approach directly. The first one
relates to the fact that the value functions is frequently merely locally
Lipschitz continuous, and the second one is that solving (2.4) suffers from
the curse of dimensionality. As was mentioned in the introduction, for this
reason many authors have tried to overcome this difficulty by using machine
learning techniques.

##  3 Methods for the construction of approximating optimal feedback laws

For enunciating precisely our goal and main results, we need first some
definitions and hypotheses, which were previously introduced in [46, 48, 45].
Throughout let $\Omega$ be an open and bounded subset of $\mathbb{R}^{d}$ and
let $\omega\subset\mathbb{R}^{d}$ be another open set with closure contained
in $\Omega$.

######  Definition 3.1.

For $u\in W^{1,\infty}_{loc}(\Omega;\mathbb{R}^{m})$, $y_{0}\in\Omega$, and
$T\in(0,\infty]$ we call $y(\cdot;y_{0},u)\in
W^{1,\infty}(0,T;\mathbb{R}^{d})$ the unique solution to (2.2) on $[0,T]$ if
it is the only function satisfying

$y^{\prime}(t)=f(y(t))+Bu(y(t)),\ y(0)=y_{0},\ \text{ for }t\in[0,T].$ (3.1)

######  Hypothesis 3.1.

Let $(T,u,\delta,y_{0})$ be a tuple with $T>0$, $u\in
W^{1,\infty}_{loc}(\Omega;\mathbb{R}^{m})$, $\delta>0$, and $y_{0}\in\Omega$.
Then the solution $y(\cdot;y_{0},u)$ exists on $[0,T]$ and
$y(t;y_{0},u)\in\Omega_{\delta}$ for $t\in[0,T)$ , where

$\Omega_{\delta}:=\\{x\in\Omega:\min_{z\in\partial\Omega}|x-z|>\delta\\}.$
(3.2)

.

######  Definition 3.2.

Let $(T,u,\delta,y_{0})$ be a tuple with $T>0$, $u\in
W^{1,\infty}_{loc}(\Omega;\mathbb{R}^{m})$, $\delta>0$, and $y_{0}\in\Omega$,
and let 3.1 hold. We define $\mathcal{V}_{T}(y_{0},u)$ as

$\mathcal{V}_{T}(y_{0},u)=\int_{0}^{T}\left(\ell(y(t;y_{0},u))+\frac{\beta}{2}|u(y(t;y_{0},u))|^{2}\right)dt.$

######  Definition 3.3.

For $v\in C^{1,1}(\Omega)$ the feedback-law $u_{v}\in
W^{1,\infty}_{loc}(\Omega;\mathbb{R}^{d})$ is defined by

$u_{v}(y)=-\frac{1}{\beta}B^{\top}\nabla v(y)\mbox{ for }y\in\Omega.$ (3.3)

We study three methods with the aim to construct an approximately optimal
feedback law. For all of them we consider a finite time horizon $T>0$, a
finite dimensional Banach space
$(\Theta,\left\lVert\cdot\right\rVert_{\Theta})$, a continuous and coercive
function $\mathcal{P}:\Theta\mapsto\mathbb{R}^{+}$, a penalty coefficient
$\alpha>0$, and a function $v:\Theta\mapsto C^{2}(\Omega)$. We will name the
tuple
$S=((\Theta,\left\lVert\cdot\right\rVert_{\Theta}),\mathcal{P},\alpha,v)$ the
setting of the learning problems.

The function $v(\theta)$ will play the role of the approximation to the value
function $V$ which will henceforth be assumed to exist at least as an element
in $L^{1}(\Omega)$. The motivation for defining the approximating problems
over finite horizons is two-fold. First, it is a natural setting for numerical
realizations, and secondly, questions related to stabilization which arise in
the context of (2.1)-(2.2) - but which are not the subject of this research -
are contained in the assumption on the well-posedness of the value functions
itself and need not be explicitly addressed for the approximating problems.

The first method, which we refer to as Averaged Feedback Learning Scheme
(AFLS), consists in minimizing the average cost of (2.1) plus a penalization
term:

######  Method 3.1.

For a given time horizon $T>0$ and a setting
$S=((\Theta,\left\lVert\cdot\right\rVert_{\Theta}),\mathcal{P},\alpha,v)$, we
solve

$\min_{\theta\in\overline{\mathcal{O}_{T}}}\frac{1}{|\omega|}\int_{\omega}\mathcal{V}_{T}(y_{0},u_{v(\theta)})dy_{0}+\alpha\mathcal{P}(\theta)$
(3.4)

where

$\mathcal{O}_{T,\Theta}:=\\{\theta\in\Theta:\mbox{ there exists
}\delta>0\mbox{ such that }u_{v_{(\theta)}}\mbox{ satifies
\lx@cref{creftypecap~refnum}{hypo:stab}}\mbox{ for all }y_{0}\in\omega\\}.$

The second method consist in finding $\theta\in\Theta$ such that it minimizes
the $L^{2}$ distance with respect to the optimal controls:

######  Method 3.2.

For $V\in C^{1}(\Omega)$, a given time horizon $T>0$, and a setting
$S=((\Theta,\left\lVert\cdot\right\rVert_{\Theta}),\mathcal{P},\alpha,v)$ we
solve

$\min_{\theta\in{\Theta}}\frac{1}{2\beta|\omega|}\int_{\omega}\int_{0}^{T}|B^{\top}(\nabla
v(\theta)(y^{*}(t;y_{0}))-\nabla
V(y^{*}(t;y_{0})))|^{2}dtdy_{0}+\alpha\mathcal{P}(\theta),$ (3.5)

where $y^{*}\in W^{1,\infty}(\omega\times(0,T);\mathbb{R}^{d})$ is such that
for almost all $y_{0}\in\omega$, the control given by
$u(t)=-\frac{1}{\beta}B^{\top}\nabla V(y^{*}(y_{0},t))$ is an optimal solution
of (2.1).

In a similar way, we consider the problem of minimizing the optimal control in
the whole domain instead of doing it along the trajectories. This leads to the
following formulation:

######  Method 3.3.

For $V\in W^{1,2}(\Omega)$ and a given setting
$S=((\Theta,\left\lVert\cdot\right\rVert_{\Theta}),\mathcal{P},\alpha,v)$ we
solve

$\min_{\theta\in\Theta}\frac{1}{2\beta|\Omega|}\int_{{\Omega}}|B^{\top}\left(\nabla
v(\theta)(y_{0})-\nabla V(y_{0})\right)|^{2}dy_{0}+\alpha\mathcal{P}(\theta).$
(3.6)

Each of the three methods is investigated for a sequence of settings
$S_{n}=((\Theta_{n},\left\lVert\cdot\right\rVert_{n}),\mathcal{P}_{n},\alpha_{n},v_{n})$
replacing $S$ in the statements of these methods, with $\Theta_{n}$ finite
dimensional. Appropriate conditions on
$S_{n}=((\Theta_{n},\left\lVert\cdot\right\rVert_{n}),\mathcal{P}_{n},\alpha_{n},v_{n})$
will be established to allow an asymptotic analysis as $n\to\infty$.

Observe that 3.1 and 3.2 when applied to a setting $S_{n}$ provide a solution
$v_{n}(\theta_{n}^{*})$ which involves initial conditions in $\omega$. But
$v_{n}(\theta_{n}^{*})$ is also defined on all of $\Omega$ and we hope that it
provides an approximation to the sought value function also outside of
$\omega$. In practice, of course, initial conditions will not be chosen in all
of $\omega$, but only at discrete sample points. These comments refer to the
topic of generalization in the machine learning context. For a specific case
of a feedback control problems, this was investigated in [48, Section 8].

######  Remark 3.1.

Methods 3.2 and 3.3 minimize the $L^{2}$ distance with respect to the optimal
feedback-law. While the first one explicitly involves information along all of
the optimal trajectories, the second one utilizes the whole domain $\Omega$.
Next we provide an interesting relation between these two formulations.
Namely, assuming that $\frac{1}{|det(D_{y_{0}}y^{*}(\cdot,\cdot))|}\in
L^{\infty}((0,T)\times\omega)$, setting
$\Gamma_{T}=\\{z\in\mathbb{R}^{d}:\exists(t,y_{0})\in[0,T]\times\omega,y^{*}(t,y_{0})=z\\}$,
and defining ${{y^{0}}}(t,z)$ for $(t,z)\in[0,T]\times\Gamma_{T}$ as the
inverse function of $y^{*}$ with respect to its second argument, we obtain by
means of the change of variables $z=y^{*}(t,y_{0})$ that

$\int_{\omega}\int_{0}^{T}|B^{\top}(\nabla v(\theta)(y^{*}(t,y_{0}))-\nabla
V(y^{*}(t,y_{0})))|^{2}dtdy_{0}=\int_{\Gamma_{T}}|B^{\top}(\nabla
v(\theta)(z)-\nabla V(z))|^{2}g(z)dz$

where

$g(z)=\int_{\\{t:\in[0,T]:\exists y_{0}\in\omega\mbox{ such that
}y^{*}(t;y_{0})=z\\}}\frac{1}{|det(D_{y_{0}}y^{*}(t,{{y^{0}}}(t,z)))|}dt.$

This reflects that the main difference between Problems 3.2 and 3.3 consists
in how they weigh the initial conditions. We also note that as consequence of
classical ODE theory the hypotheses on ${det(D_{y_{0}}y^{*})}$ expressed above
hold if $V\in C^{2}(\Omega)$ .

To address the problem of approximating $v$ in Methods 3.1-3.3 by sequences of
finite dimensional settings
$S_{n}=((\Theta_{n},\left\lVert\cdot\right\rVert_{n}),\mathcal{P}_{n},\alpha_{n},v_{n})$
we require additional hypotheses which are explained next. This will need some
patience on the side of the reader. Let us consider sequences of time horizons
$T_{n}>0$, of finite dimensional Bachach spaces
$\\{(\Theta_{n},\left\lVert\cdot\right\rVert_{n}\\}_{n\in\mathbb{N}}$, penalty
operators $\mathcal{P}_{n}:\Theta_{n}\mapsto\mathbb{R}$, functions $v_{n}\in
C^{1}(\Theta_{n};C^{2}(\Omega))$, and penalty parameters $\alpha_{n}>0$. Let
$\theta_{n}^{*}\in\Theta_{n}$ be obtained by solving any of the three problems
introduced above. Our goal then consists in finding assumptions which imply
that

$\lim_{n\to\infty}\left\lVert\mathcal{V}_{T_{n}}(\cdot,u_{v_{n}(\theta_{n}^{*})})-V\right\rVert_{L^{p}(\omega)}=0,$
(3.7)

for $p\in[1,\infty]$. For this, we will use the following hypotheses:

######  Hypothesis 3.2.

We have $V\in L^{1}(\Omega)$ and there exist an increasing sequence
$\hat{T}_{n}$ of time horizons, a constant
$\delta\in(0,dist(\omega,\partial\Omega))$, and a sequence $V_{n}\in
C^{2}(\Omega)$ such that $\lim_{n\to\infty}\hat{T}_{n}=\infty$, for all
$y_{0}\in\omega$ the tuple $(\hat{T}_{n},u_{V_{n}},\delta,y_{0})$ satisfies
3.1 with $T=\hat{T}_{n}$, and

$\limsup_{n\to\infty}\int_{\omega}(\mathcal{V}_{\hat{T}_{n}}(y_{0},u_{V_{n}})-V(y_{0}))dy_{0}=0.$
(3.8)

While the previous hypothesis serves the purpose of guaranteeing that $V$ can
be approximated by sufficiently smooth functions, the following one will
ensure the possibility to approximate these smooth functions by finitely
parameterized functions.

######  Hypothesis 3.3.

Let
$\\{(\Theta_{n},\left\lVert\cdot\right\rVert_{\Theta_{n}}\\}_{n\in\mathbb{N}}$
be a sequence of finite dimensional Banach spaces and $v_{n}:\Theta_{n}\mapsto
C_{loc}^{2}(\Omega)$ be continuous functions. For each $g\in C^{1}(\Omega)$
there exist $\theta_{n}$ and a sequence of penalty coefficients $\alpha_{n}>0$
with $\lim_{n\to\infty}\alpha_{n}=0$ such that

$\lim_{n\to\infty}v_{n}(\theta_{n})=g\mbox{ in }C_{loc}^{1}(\Omega)\mbox{ and
}\lim_{n\to\infty}\alpha_{n}\mathcal{P}_{n}(\theta_{n})=0.$ (3.9)

Hypotheses 3.2 and 3.3, and a proper choice of the time horizons $T_{n}$ will
imply the convergence of 3.1 in the sense of (3.7). Regarding the feasibility
of 3.2, it was proved in [46] that this hypothesis holds if $V$ is
$\alpha-$Hölder continuous with $\alpha>\frac{1}{2}$ and that there exist
$w\in C^{1}(\Omega)$ and $\delta>0$ such that 3.5 is satisfied with $\phi=V$.

A more restrictive version of 3.3, requiring higher regularity of the value
function, will allow us to prove that Methods 3.2 and 3.3 converge in the same
manner as specified for 3.1 above.

######  Hypothesis 3.4.

There exists a Banach space $(\Theta,\left\lVert\cdot\right\rVert_{\Theta})$
such that $V\in\Theta\cap C^{2}(\Omega)$, $\Theta$ is compactly embedded in
$C^{1}(\overline{\Omega})$ and there exists $\delta>0$ such that for all
$y_{0}\in\omega$ the tuple $(T,u_{V},\delta,y_{0})$ satisfies 3.1 with
$T=\infty$. Further, there exist sequences
$\\{(\Theta_{n},\left\lVert\cdot\right\rVert_{n})\\}_{n\in\mathbb{N}}$ of
finite dimensional Banach spaces, $v_{n}:\Theta_{n}\mapsto\Theta$ of
continuous functions, and $P_{n}:\Theta_{n}\mapsto\Theta$, such that for all
$g\in\Theta$ there exists a sequence $\theta_{n}\in\Theta_{n}$ satisfying

$\lim_{n\to\infty}\left\lVert v_{n}(\theta_{n})-g\right\rVert_{\Theta}=0,\
\sup_{n\in\mathbb{N}}\mathcal{P}_{n}(\theta_{n})<\infty.$

Further there exist constants $C>0$ and $\sigma>0$ such that it holds

$\left\lVert v_{n}(\theta)\right\rVert_{\Theta}\leqslant
C\mathcal{P}_{n}(\theta)^{\sigma},\text{ for all }n\in\mathbb{N}\text{ and
}\theta\in\Theta_{n}.$

######  Remark 3.2.

In order to illustrate the notation and give a situation where 3.4 holds, we
consider the case of polynomials as ansatz functions. For the sake of
simplicity we take $\Omega=(-1,1)^{d}$. For $k\in\mathbb{N}$ consider
$H^{k}([-1,1])$ with the inner product
$\langle\phi,\psi\rangle=\frac{1}{2}\sum_{i=0}^{k}\int_{-1}^{1}\phi^{(i)}(x)\psi^{(i)}(x)dx$
where $\phi^{i}$ stands for the $i$-th derivative of $\phi$. Let
$\\{\phi^{k}_{i}\\}_{i=1}^{\infty}$ be a orthonormal basis of $H^{k}((-1,1))$
obtained by applying the Gram-Schmidt process to the set of monomials
$\\{1,x,x^{2},\ldots\\}$. We set
$\Theta=H_{mix}^{k}((-1,1)^{d}):=\otimes_{i=1}^{d}H^{k}((-1,1))$, i.e,
$\Theta$ is the tonsorial product of $d$ copies of $H^{k}((-1,1))$. The space
$H_{mix}^{k}((-1,1)^{d})$ is Hilbert with the tensorial inner product. For a
multi-index $\alpha\in\mathbb{N}^{d}$ we denote
$\phi_{\alpha}(y):=\prod_{i=1}^{d}\phi_{i}(y_{i})^{\alpha_{i}}$. Hence the set
$\\{\phi_{\alpha}\\}_{\alpha\in\mathbb{N}^{d}}$ is a basis for $\Theta$. For
$n\in\mathbb{N}$ we define $\Lambda_{n}\subset\mathbb{N}^{d}$ as the follows:

$\Lambda_{n}=\\{\alpha=(\alpha_{1},\ldots,\alpha_{d})\in\mathbb{N}^{d}:\max_{i=1,\ldots,d}\alpha_{i}\leqslant
n\\}.$

Let us consider an enumeration of $\Lambda_{n}$ given by
$\\{\alpha^{i}\\}_{i=1}^{|\Lambda_{n}|}$. We set
$\Theta_{n}=\mathbb{R}^{|\Lambda_{n}|}$,
$v_{n}(\theta)=\sum_{i=1}^{|\Lambda_{n}|}\theta_{i}\phi_{\alpha^{i}}$,
$\left\lVert\theta\right\rVert_{n}=|\theta|$, and
$\mathcal{P}_{n}(\theta)=\frac{1}{2}|\theta|^{2}$. With these choices 3.4 is
met with $\sigma=\frac{1}{2}$ and $k\geqslant 3$, since
$H_{min}^{k}((-1,1)^{d})$ is continuously embedded in
$W^{2,\infty}((-1,1)^{d})$ (see [49]).

For high-dimensional problems (i.e. large $d$) the hyperbolic-cross index set
is a more efficient choice than $\Lambda_{n}$. This set of multi-indexes is
given by

$\Gamma_{n}=\\{\alpha=(\alpha_{1},\ldots,\alpha_{d})\in
N^{d}:\prod_{i=1}^{d}(\alpha_{i}+1)\leqslant n\\}.$

For $\Omega=\otimes_{i=1}^{d}(a_{i},b_{i})$ with $a_{i}<b_{i}$ for
$i\in\\{1,\ldots,d\\}$ a change of variables is needed, that is, for
$x\in\Omega$ we define
$\tilde{\phi}_{\alpha}(x)=\prod_{i=1}^{d}\phi_{\alpha_{i}}\left(\frac{2}{(b_{i}-a_{i})}\left(x_{i}-\frac{a_{i}+b_{i}}{2}\right)\right)$
and for $\theta\in\mathbb{R}^{|\Lambda_{n}|}$ we set
$\tilde{v}_{n}(\theta)=\sum_{i=1}^{n}\theta_{i}\tilde{\phi}_{\alpha^{i}}$.
With reference to 3.4 we choose $\Theta_{n}=\mathbb{R}^{|\Lambda_{n}|}$ and
$\Theta=H^{k}_{mix}((-1,1)^{d})$ as above. Again the embedding of $\Theta$
into $C^{1}(\bar{\Omega})$ is compact for $k\geqslant 3$.

In a more general perspective, 3.4 also encsompasses other interesting cases
such as neural networks with $\Theta$ being the respective Barron space
provided the activation function is smooth, or reproducing kernel Hilbert
spaces with a $C^{2}$ kernel and $\Theta$ as the corresponding native space.
In such cases $\Theta$ can then be chosen as the space of limits of functions
of the form $v_{n}(\theta_{n})$ for some $\theta_{n}\in\Theta_{n}$, endowed
with the proper norm.

Let us now return to 3.1. We next state an alternative hypothesis involving a
Ljapunov function $w$ and in the subsequent remark we discuss that this
hypothesis implies 3.1. Recall that $\omega$ is an open set satisfying
$\omega\Subset\Omega$.

######  Hypothesis 3.5.

There exist $\delta>0$, and $w\in C^{1}(\Omega)$ such that $\omega$ is
strictly contained in

$\omega_{\delta}:=\\{y\in\Omega:w(y)<\sup_{y_{0}\in\omega}w(y_{0})+\delta\\},$

$\overline{\omega_{\delta}}$ is strictly contained in $\Omega$, and
$\partial\omega_{\delta}$ is of class $C^{1}$. Moreover $\phi\in C(\Omega)$ is
a viscosity super solution of

$-\nabla w(y)^{\top}(f(y)+B(y)u_{\phi}(y))=0\mbox{ in }\omega_{\delta},$

i.e. for every $\bar{y}\in\omega_{\delta}$ and every $h\in
C_{loc}^{1}(\Omega)$ such that $\phi-h$ attains a local minimum at $\bar{y}$
the inequality

$\nabla w(\bar{y})^{\top}(f(\bar{y})+B(\bar{y})u_{h}(\bar{y}))\leqslant 0$

holds.

######  Remark 3.3.

Let us consider $\delta>0$, $w\in C^{1}$, $\phi\in C(\Omega)$ as in 3.5, and
$T>0$. We note that under these conditions we have that for all
$y_{0}\in\omega$ the tuple $(T,u_{\phi},\sigma,y_{0})$ satisfies 3.1, since
$\overline{\omega_{\delta}}$ is strictly included in $\Omega$, where $\sigma$
is the distance of $\omega_{\delta}$ to $\partial\Omega$, so that
$\omega_{\delta}\subset\Omega_{\sigma}$.

We also recall from Lemma 7.1 in [46] that if $\phi\in C(\Omega)$ satisfies
3.5, it is possible to choose $\sigma>0$ and to construct a sequence
$(\phi_{n})_{n\in\mathbb{N}}$ in $C^{2}(\Omega)$ such that for all
$n\in\mathbb{N}$ and for all $y_{0}\in\omega$ the tuple
$(T,u_{\phi_{n}},\sigma,y_{0})$ satisfies 3.1 and 3.2.

######  Remark 3.4.

Note that in contrast to 3.1, knowledge of the value function is required for
3.2 and 3.3. In practice this can be realized by approximating the integrals
by finite sums or by Monte Carlo methods with values for $\nabla V$ at sample
points obtained by open loop control solves. This approach was investigated in
[22] for instance. Recall the $\nabla V$ is related to the adjoint state. It
is obtained as a by-product in every numerical approach for open loop optimal
control. But, as already pointed out above, approximating the feedback over
the smaller set $\omega$ by globally defined functions also serves the purpose
to utilize the approximating feedback laws outside of $\omega$. This
generalization step, of course, requires additional analysis and will depend
on the underlying dynamics characterized by $f$.

##  4 Convergence analysis

This section is devoted to the convergence analysis for Methods 3.1 \- 3.3. We
first state the two main results and then provide their proofs below. For 3.1
we have the following convergence property.

######  Theorem 4.1.

Let
$S_{n}=((\Theta_{n},\left\lVert\cdot\right\rVert_{n}),\mathcal{P}_{n},\alpha_{n},v_{n})$
be a sequence of settings such that 3.3 holds and assume 3.2 holds as well.
Then there exists a sequence of times $T_{n}>0$ with
$\lim_{n\to\infty}T_{n}=\infty$, and a sub-sequence of settings $S_{n}$
(denoted by the same index) for which 3.1 with $S=S_{n}$ has at least one
solution. Further, let $\theta^{*}_{n}\in\overline{\mathcal{O}_{n}}$ be a
sequence of solutions to 3.1 with $S=S_{n}$. Then we have

$\lim_{n\to\infty}\mathcal{V}_{T_{n}}(\cdot,u_{v_{n}(\theta_{n})})=V(\cdot)\mbox{
in }L^{1}(\omega).$

In order to prove the convergence for 3.3, we shall need to assume 3.4 instead
of Hypotheses 3.2 and 3.3. In this case, the convergence holds in
$L^{\infty}(\omega)$. The regularity assumption on $V$ expressed in Hypothesis
3.4, however, is more restrictive than the hypotheses for Theorem 4.1 where
the assumptions on $V$ are only those involved in Hypothesis 3.2 . We shall
see that the convergence for 3.3 also implies the convergence for 3.2, by
using a change of variables as was already mentioned in 3.1.

######  Theorem 4.2.

Let $V$, $(\Theta,\left\lVert\cdot\right\rVert_{\Theta})$, and a sequence of
settings
$S_{n}=((\Theta_{n},\left\lVert\cdot\right\rVert_{n}),\mathcal{P}_{n},\alpha_{n},v_{n})$
be such that 3.4 holds, and let $\lim_{n\to\infty}\alpha_{n}=0$. Then there
exists a sequence of times $T_{n}>0$ for which every sequence
$\theta^{*}_{n}\in\overline{\mathcal{O}_{T_{n}}}$ of solutions of either 3.2
or 3.3 with $S=S_{n}$ and $T=T_{n}$ satisfies that
$\theta_{n}^{*}\in\mathcal{O}_{T_{n}}$ and

$\lim_{n\to\infty}\mathcal{V}_{T_{n}}(\cdot,u_{v_{n}(\theta_{n})})=V\mbox{ in
}L^{\infty}(\omega).$

To verify these two results we start by proving the equivalence between (3.8),
and the following two properties

$\lim_{n\to\infty}\mathcal{V}_{\hat{T}_{n}}(\cdot,u_{V_{n}})=V\mbox{ in
}L^{1}(\omega),$ (4.1)

and

$\lim_{n\to\infty}\int_{\omega}(\mathcal{V}_{\hat{T}_{n}}(y_{0},u_{V_{n}})-V(y_{0}))^{+}dy_{0}=0.$
(4.2)

######  Lemma 4.1.

Let $T_{n}>0$ be an increasing sequence with $\lim_{n\to\infty}T_{n}=\infty$,
let $\delta$ be a positive constant, and let $V_{n}\in
C^{2}(\overline{\Omega})$ be as sequence, such that the tuple
$(T_{n},u_{V_{n}},\delta,y_{0})$ satisfies 3.1 for all $y_{0}\in\omega$ and
$n\in\mathbb{N}$. Then for all $y_{0}\in\omega$ the following identity holds

$\lim_{n\to\infty}(V(y_{0})-\mathcal{V}_{T_{n}}(y_{0},u_{V_{n}}))^{+}=0.$
(4.3)

###### Proof.

By contradiction, let assume there exists $y_{0}$ such that (4.3) does not
hold. This implies that there exist $\kappa>0$ and a sub-sequence still
denoted by $V_{n}$ such that

$(V(y_{0})-\mathcal{V}_{T_{n}}(y_{0},u_{V_{n}}))^{+}>\kappa$ (4.4)

for all $n\in\mathbb{N}$. Let us denote $y_{n}=y(t;y_{0},u_{V_{n}})$ and
$u_{n}=u_{V_{n}}\circ y_{n}$. By 3.1 and the boundedness of $\Omega$ we have
that $\\{\|y_{n}\|_{L^{\infty}((0,T_{n});\mathbb{R}^{d})}\\}_{n=1}^{\infty}$
is bounded. By (4.4) and the definition of $u_{n}$ we have that for all
$n\in\mathbb{N}$

$\int_{0}^{T_{n}}\left(\ell(y_{n}(t))+\frac{\beta}{2}|u_{n}(t)|^{2}\right)dt\leqslant
V(y_{0})-\kappa.$ (4.5)

This implies that
$\\{\|u_{n}\|_{L^{2}((0,T_{n});\mathbb{R}^{m})}\\}_{n=1}^{\infty}$ is bounded.
This, together with the fact that $y_{n}$ is a solution of

$\frac{d}{dt}y_{n}(t)=f(y_{n}(t))+Bu_{n}(t),\ t\in(0,T_{n}),\ y_{n}(0)=y_{0},$
(4.6)

that $f$ is Lipschitz continuous on bounded sets, and
$\\{y_{n}\\}_{n=1}^{\infty}$ is bounded in
$L^{\infty}((0,T_{n});\mathbb{R}^{d})$, we obtain that
$\\{y_{n}\\}_{n=1}^{\infty}$ is uniformly bounded in
$H^{1}((0,T_{n});\mathbb{R}^{d})$. Therefore, by a diagonalization argument
there exists $y^{*}\in H_{loc}^{1}((0,\infty);\mathbb{R}^{d})$ and $u\in
L_{loc}^{2}((0,\infty);\mathbb{R}^{m})$ such that by passing to a sub-
sequence, if needed, we have

$u_{n}\rightharpoonup u^{*}\mbox{ in }L_{loc}^{2}((0,\infty);\mathbb{R}^{m}),\
y_{n}\rightharpoonup y^{*}\mbox{ in
}H_{loc}^{1}((0,\infty);\mathbb{R}^{d}),\mbox{ and }y_{n}\to y^{*}\mbox{ in
}C_{loc}([0,\infty);\mathbb{R}^{d}),$ (4.7)

as $n$ tends to infinity. Using these facts together with (4.6) we have that

$\frac{d}{dt}y^{*}(t)=f(y^{*}(t))+Bu^{*}(t),\ t\in(0,\infty),\
y^{*}(0)=y_{0}.$ (4.8)

Next, let $T>0$ be a fixed time horizon and choose $n_{0}\in\mathbb{N}$ such
that for all $n>n_{0}$ we have $T<T_{n}$. Hence by (4.5) we get that for all
$n>n_{0}$

$\int_{0}^{T}\left(\ell(y_{n}(t))+\frac{\beta}{2}|u_{n}(t)|^{2}\right)dt\leqslant
V(y_{0})-\kappa.$ (4.9)

Using (4.7), the continuity of $\ell$ and the weak lower semi-continuity of
the norm of $L^{2}((0,T);\mathbb{R}^{m})$, letting $n$ tend to infinity in
(4.9) we obtain

$\int_{0}^{T}\left(\ell(y^{*}(t))+\frac{\beta}{2}|u^{*}(t)|^{2}\right)dt\leqslant
V(y_{0})-\kappa.$ (4.10)

Since this holds for any arbitrary $T$ we get that $u^{*}\in
L^{2}((0,\infty);\mathbb{R}^{d})$ and

$\int_{0}^{\infty}\left(\ell(y^{*}(t))+\frac{\beta}{2}|u^{*}(t)|^{2}\right)dt\leqslant
V(y_{0})-\kappa,$

which contradicts the fact that $V$ is the value function. Consequently, we
have that for all $y_{0}\in\omega$

$\lim_{n\to\infty}(V(y_{0})-\mathcal{V}_{T_{n}}(y_{0},u_{V_{n}}))^{+}=0.$

∎

######  Corollary 4.1.

Under the hypotheses of the previous lemma (3.8), (4.1) and (4.2) are
equivalent.

###### Proof.

Since $(V(y_{0})-\mathcal{V}_{T_{n}}(y_{0},u_{V_{n}}))^{+}\leqslant V(y_{0})$
we can apply the dominated convergence theorem together with (4.3) to get that

$\lim_{n\to\infty}\int_{\omega}(V(y_{0})-\mathcal{V}_{T_{n}}(y_{0},u_{V_{n}}))^{+}dy_{0}=0.$
(4.11)

Rewriting the $L^{1}(\omega)$ norm of
$V(y_{0})-\mathcal{V}_{T_{n}}(y_{0},u_{V_{n}})$ we find that

$\int_{\omega}|V(y_{0})-\mathcal{V}_{T_{n}}(y_{0},u_{V_{n}})|dy_{0}=\int_{\omega}(\mathcal{V}_{T_{n}}(y_{0},u_{V_{n}}))-V(y_{0}))^{+}dy_{0}+\int_{\omega}(V(y_{0})-\mathcal{V}_{T_{n}}(y_{0},u_{V_{n}}))^{+}dy_{0}.$

Since the convergence to 0 of the last term in the previous inequality is
always granted under the hypotheses, we get the equivalence of (4.1) and
(4.2). Further, condition (4.1) clearly implies (3.8). To conclude we only
need to prove that (3.8) implies either (4.1) or (4.2). For this purpose we
note that

$\begin{array}[]{ll}0&\leqslant\liminf\int_{\omega}(\mathcal{V}_{T_{n}}(y_{0},u_{V_{n}})-V(y_{0}))^{+}dy_{0}\\\\[6.45831pt]
&\leqslant\lim\int_{\omega}(V(y_{0})-\mathcal{V}_{T_{n}}(y_{0},u_{V_{n}}))^{+}dy_{0}+\limsup\int_{\omega}(\mathcal{V}_{T_{n}}(y_{0},u_{V_{n}})-V(y_{0}))dy_{0}.\end{array}$

Hence, if (3.8) holds, we can combine it with (4.11) to get (4.2).

∎

Before continuing we need to establish the following technical result
concerning the stability of the solution to (3.1) if $u$ is of the form (3.3)
and perturbations occur in $v$. Its proof uses Lemma A1.

######  Corollary 4.2.

Let $T>0$ be a time horizon, $y_{0}\in\omega$, $g_{1}\in C^{2}(\Omega)$,
$g_{2}\in C^{2}(\Omega)$ and $\delta>0$, be such that
$(T,u_{g_{2}},\delta,y_{0})$ satisfies 3.1. Assume that

$\frac{|B|\left\lVert
u_{g_{1}}-u_{g_{2}}\right\rVert_{L^{\infty}(\Omega_{\frac{\delta}{4}};\mathbb{R}^{m})}}{a}\left(e^{Ta}-1\right)\leqslant\frac{\delta}{2},$
(4.12)

where

$a=\left\lVert
Df\right\rVert_{L^{\infty}(\Omega_{\frac{\delta}{4}};\mathbb{R}^{d\times
d})}+\frac{|B|^{2}}{\beta}\left\lVert\nabla^{2}g_{2}\right\rVert_{L^{\infty}(\Omega_{\frac{\delta}{4}};\mathbb{R}^{d\times
d})}.$ (4.13)

Then $(T,u_{g_{1}},\frac{\delta}{2},y_{0})$ satisfies 3.1.

###### Proof.

For the sake abbreviation we set $y_{1}(t)=y(t;y_{0},u_{g_{1}})$ and
$y_{2}(t)=y(t;y_{0},u_{g_{2}})$ for all $t\geqslant 0$ such that these
expressions are well defined. Let $\hat{T}>0$ be the escape time of $y_{1}$
from $\Omega_{\frac{\delta}{2}}$, i.e.
$\hat{T}=\sup\\{\tilde{T}>0:y_{1}(t)\in\Omega_{\frac{\delta}{2}}\mbox{ for all
}t\in[0,\tilde{T}]\\}$. We need to prove that $\hat{T}\geqslant T$. By
contradiction, assume that $\hat{T}<T$. Then there exists
$\tilde{T}\in(\hat{T},T)$ such that for all $t\in(\hat{T},\tilde{T}]$ we have
$y_{1}(t)\in\Omega_{\frac{\delta}{4}}\setminus\Omega_{\frac{\delta}{2}}$.
Combining (A.1) ( with $\delta$ replaced by $\frac{\delta}{4}$) and (4.12) we
obtain for all $t\in(\hat{T},\tilde{T}]$ that

$|y_{1}(t)-y_{2}(t)|\leqslant\frac{\delta}{2}.$

Here we use that if $(T,u_{g_{2}},\delta,y_{0})$ satisfies 3.1 then it is also
satisfied with $(T,u_{g_{2}},\frac{\delta}{4},y_{0})$. Further again using
that $(T,u_{g_{2}},\delta,y_{0})$ satisfies 3.1 and that $\tilde{T}<T$ we get
for all $t\in(\hat{T},\tilde{T}]$ and for all $z\in\partial\Omega$

$\delta\leqslant|z-y_{2}(t)|<|z-y_{1}(t)|+|y_{1}(t)-y_{2}(t)|\leqslant\frac{\delta}{2}+|z-y_{1}(t)|.$

This implies that for $t\in(\hat{T},\tilde{T}]$ we have
$|z-y_{1}(t)|>\frac{\delta}{2}$ for all $z\in\partial\Omega$ which implies
that $y_{1}(t)\in\Omega_{\frac{\delta}{2}}$ for $t\in(\hat{T},\tilde{T}]$.
This contradicts the definition of $\hat{T}$ and therefore $\hat{T}\geqslant
T$. ∎

###### Proof of Theorem 4.1.

For proving the existence of the sub-sequence of $S_{n}$ with
$\lim_{n\to\infty}T_{n}=\infty$ such that 3.1 has a solution we proceed
similarly to the proof of Lemma 3 in [48]. Let $T_{n}$, $V_{n}\in
C^{2}(\Omega)$ and $\delta$ be the quantities appearing in 3.2. Next we use
3.3 which determines $(\Theta_{k},\|\cdot\|_{\Theta_{k}})$ and $v_{k}$, and
choose with $g=V_{n}$, for each $n\in\mathbb{N}$ a sequence
$\theta_{k,n}\in\Theta_{k}$ such that
$\lim_{k\to\infty}v_{k}(\theta_{k,n})=V_{n}\mbox{ in }C_{loc}^{1}(\Omega)$ and
$\lim_{k\to\infty}\alpha_{k}\mathcal{P}_{k}(\theta_{k,n})=0$.

Let $y_{0}\in\omega$ be fixed and set $u_{n}(x)=u_{V_{n}}(x)$ and
$u_{k,n}(x)=u_{v_{k}(\theta_{k,n})}(x)$ for all $x\in\Omega$ and
$n,k\in\mathbb{N}$. Let
$L_{n}=\left\lVert\nabla^{2}V_{n}\right\rVert_{L^{\infty}\big{(}\Omega_{\frac{\delta}{4}};\mathbb{R}^{d\times
d}\big{)}}$. Since for each $n\in\mathbb{N}$, $v_{k}(\theta_{k,n})$ converges
to $V_{n}$ in $C_{loc}^{1}(\Omega)$, there exists $k(n)$ such that for all
$k>k(n)$ we have

$\left(e^{(L_{f}+\frac{|B|^{2}}{\beta}L_{n})T_{n}}-1\right)\frac{|B|\left\lVert
u^{k,n}-u^{n}\right\rVert_{L^{\infty}(\Omega_{\frac{\delta}{4}};\mathbb{R}^{m})}}{(L_{f}+\frac{|B|^{2}}{\beta}L_{n})}\leqslant\frac{\delta}{2}\mbox{
and }\left\lVert
u^{k,n}-u^{n}\right\rVert_{L^{\infty}(\Omega_{\frac{\delta}{4}};\mathbb{R}^{m})}\leqslant
1.$ (4.14)

In view of (4.14) and the fact that $(T_{n},u_{n},\delta,y_{0})$ satisfies 3.1
we can apply 4.2 with $g_{1}=v_{k}(\theta_{k,n})$, $g_{2}=V_{n}$ and $T=T_{n}$
to obtain that for all $y_{0}\in\omega$ and $k\geqslant k(n)$ that the tuple
$(T_{n},u^{k,n},\frac{\delta}{2},y_{0})$ satisfies 3.1. Consequently, we have
that $\theta_{k,n}\in\mathcal{O}_{T_{n},\Theta_{k}}$ for $k>k(n)$. Combining
this, the fact that $\mathcal{P}_{k}$ is coercive and continuous, and that
$\mathcal{V}_{T_{n}}$ is continuous, we get by the direct method that 3.1 with
$S=S_{k}$ and $T=T_{n}$ for $k>k(n)$ has an optimal solution.

For the second part we use (A.2). With this in mind we set

$K_{n}=|B|\left(\left\lVert\nabla\ell\right\rVert_{L^{\infty}(\Omega;\mathbb{R}^{d})}\left(e^{(L_{f}+\frac{|B|^{2}}{\beta}L_{n})T_{n}}-1\right)\frac{1}{(L_{f}+\frac{|B|^{2}}{\beta}L_{n})^{2}}+\frac{T_{n}\beta}{2}(\left\lVert
u^{n}\right\rVert_{L^{\infty}(\Omega_{\frac{\delta}{2}})}+1)\right).$

Using (A.2) we obtain for all $n\in\mathbb{N}$ and $k>k(n)$

$\int_{\omega}\mathcal{V}_{T_{n}}(y_{0},v_{k}(\theta_{k,n}))dy_{0}\leqslant\int_{\omega}\mathcal{V}_{T_{n}}(y_{0},V_{n})dy_{0}+|\omega|K_{n}\left\lVert
u^{k,n}-u^{n}\right\rVert_{L^{\infty}(\Omega_{\frac{\delta}{4}};\mathbb{R}^{m})}$
(4.15)

By 3.3 we can choose $k_{n}\geqslant k(n)$ with
$\lim_{n\to\infty}k_{n}=\infty$, such that

$\lim_{n\to\infty}|\omega|K_{n}\left\lVert
u^{k_{n},n}-u^{n}\right\rVert_{L^{\infty}(\Omega_{\frac{\delta}{4}};\mathbb{R}^{m})}+\alpha_{k_{n}}\mathcal{P}_{k_{n}}(\theta_{k_{n},n})=0.$
(4.16)

Let $\theta_{n}^{*}\in\Theta_{k_{n}}$ be an optimal solution of 3.1 for
$T=T_{n}$ and $S=S_{k_{n}}$. Then by (4.15)

$\begin{array}[]{l}\displaystyle\int_{\omega}\mathcal{V}_{T_{n}}(y_{0},v_{k_{n}}(\theta^{*}_{k_{n}}))dy_{0}+\alpha_{k_{n}}\mathcal{P}_{k_{n}}(\theta_{k_{n}})\leqslant\int_{\omega}\mathcal{V}_{T_{n}}(y_{0},v_{k_{n}}(\theta_{k_{n},n}))dy_{0}+\alpha_{k_{n}}\mathcal{P}_{k_{n}}(\theta_{k_{n},n})\\\
\displaystyle\leqslant\int_{\omega}\mathcal{V}_{T_{n}}(y_{0},u_{V_{n}})dy_{0}+\alpha_{k_{n}}\mathcal{P}_{k_{n}}(\theta_{k_{n},n})+|\omega|K_{n}\left\lVert
u^{k,n}-u^{n}\right\rVert_{L^{\infty}(\Omega_{\frac{\delta}{4}};\mathbb{R}^{m})}.\end{array}$

Combining the above inequality, (4.16) , and 3.2, we arrive at

$\limsup_{n\to\infty}\int_{\omega}\mathcal{V}_{T_{n}}(y_{0},v_{k_{n}}(\theta^{*}_{k_{n}}))dy_{0}\leqslant\int_{\omega}V(y_{0})dy_{0}.$

Therefore, by 4.1 we conclude the proof. ∎

###### Proof of Theorem 4.2.

Let us start by considering $\theta^{*}_{n}$ as the solution of 3.3. From the
optimality of $\theta^{*}_{n}$ we obtain that

$\frac{1}{2\beta}\left\lVert B^{\top}(\nabla v_{n}(\theta^{*}_{n})-\nabla
V)\right\rVert_{L^{2}(\Omega)}^{2}+\alpha_{n}\mathcal{P}_{n}(\theta_{n}^{*})\leqslant\frac{1}{2\beta}\left\lVert
B^{\top}(\nabla v_{n}(\theta_{n})-\nabla
V)\right\rVert_{L^{2}(\Omega)}^{2}+\alpha_{n}\mathcal{P}_{n}(\theta_{n})$
(4.17)

with $\theta_{n}$ as in 3.4. By 3.4 we know that the right hand side of this
inequality converges to $0$ as $n$ goes to infinity and consequently we obtain
that

$\lim_{n\to\infty}B^{\top}\nabla v_{n}(\theta^{*}_{n})=B^{\top}V\mbox{ in
}L^{2}(\Omega).$

Moreover, again by 3.4 and the optimality of $\theta_{n}^{*}$ we obtain that
$\left\lVert v_{n}(\theta_{n}^{*})\right\rVert_{\Theta}$ is uniformly bounded
for $n\in\mathbb{N}$. This together with the fact that $\Theta$ is compactly
embedded in $C^{1}(\overline{\Omega})$ implies that $B^{\top}\nabla
v_{n}(\theta_{n}^{*})$ converges to $B^{\top}\nabla V$ in
$C^{1}(\overline{\Omega})$. Let $T_{n}>0$ be a sequence such that
$\lim_{n\to\infty}e^{T_{n}a}\left\lVert B^{\top}(\nabla
v_{n}(\theta_{n}^{*})-\nabla
V)\right\rVert_{L^{\infty}(\Omega_{\frac{\delta}{2}};\mathbb{R}^{m})}=0$,
where $a$ is given by (4.13) with $g_{2}=V$. Using 4.2 we know that for all
$y_{0}\in\omega$ the tuple
$(T_{n},u_{v_{n}(\theta_{n}^{*})},\frac{\delta}{2},y_{0})$ satisfies 3.1 with
$\delta$ replaced by $\frac{\delta}{2}$ for $n$ large enough. With Lemma A1 we
conclude that

$\lim_{n\to\infty}\mathcal{V}_{T_{n}}(\cdot,u_{v_{n}(\theta_{n}^{*})})=V\mbox{
in }L^{\infty}(\omega)$

and thus the claim for 3.3 is verified.

Next we turn now our attention to 3.2. Let $\theta_{n}^{*}$ be an optimal
solution for 3.2 with $S=S_{n}$ and $T=T_{n}$ as before. Since
$D_{y_{0}}y^{*}(t,y_{0})$ is invertible for all
$(t,y_{0})\in[0,T]\times\overline{\omega}$, by using the change of variables
$z=y^{*}(t;y_{0})$ as in 3.1 we obtain

$\int_{\omega}\int_{0}^{T}|B^{\top}(\nabla v(y^{*}(t;y_{0}))-\nabla
V(y^{*}(t;y_{0})))|^{2}dtdy_{0}=\int_{\Gamma_{T}}|B^{\top}(\nabla v(z)-\nabla
V(z))|^{2}g(z)dz$

Moreover, we know that $g$ is strictly positive in $\Gamma_{T}$ and therefore
it is equivalent to the Lebesgue measure in $\Gamma_{T}$. From this we can
continue with the proof as in the case of 3.3 to obtain the desired
convergence. ∎

######  Remark 4.1.

Here we informally provide connections between 3.1 and 3.2 which are obtained
through their optimality conditions. Let us consider a time horizon
$T\in(0,\infty)$ and a setting
$S=((\Theta,\left\lVert\cdot\right\rVert_{\Theta}),\mathcal{P},\alpha,v)$. The
differentiability of the function $\mathcal{V}_{T}(y_{0},u)$ for
$y_{0}\in\Omega$ and $u\in C^{2}(\Omega)$ with respect to both arguments comes
from [46] (see also A). Using this and (A.8) we get the optimality condition
of 3.1 as

$\frac{1}{\beta|\omega|}\\!\int_{\omega}\int_{0}^{T}\\!\\!\big{(}\nabla_{y}v(\theta)(y(t;y_{0},u_{v(\theta)}))-\nabla\\!_{y}\\!\mathcal{V}_{T}(y(t;y_{0},u_{v(\theta)}),u_{v(\theta)}\big{)}^{\top}\\!\\!BB^{\top}\nabla\psi(y(t;y_{0},u_{v(\theta)}))dtdy_{0}+\alpha
D\mathcal{P}(\theta)\cdot\vartheta=0,$ (4.18)

where $\psi=\left(D_{\theta}v\cdot\vartheta\right)$, and $\vartheta\in\Theta$
is arbitrary. We note that (4.18) also corresponds to the optimality
conditions of the regression problem given by

$\min_{\vartheta\in\Theta}\frac{1}{2\beta|\omega|}\int_{\omega}\int_{0}^{T}|B^{\top}(\nabla_{y}v(\vartheta)(y(t;y_{0},u_{v(\theta)}))-\nabla\\!_{y}\mathcal{V}_{T}(y(t;y_{0},u_{v(\theta)}),u_{v(\theta)}))|^{2}dtdy_{0}+\alpha\mathcal{P}(\vartheta)$
(4.19)

It is noteworthy to observe that this problem is similar to the problem
arising in 3.2, replacing $V$ by $\mathcal{V}_{T}$ and $y^{*}$ by
$y(\cdot;\cdot,u_{v(\vartheta)})$. This implies that if $\theta^{*}\in\Theta$
is an optimal solution of 3.2 and the trajectories
$y(\cdot;y_{0},u_{v(\theta^{*})})$ are close to $y^{*}(\cdot;\cdot)$ in
$W^{1,\infty}((0,T)\times\omega;\mathbb{R}^{d})$, then $v(\theta^{*})$
satisfies approximately (4.18). This holds in the opposite direction as well,
that is, if $\theta^{*}$ is an optimal solution of 3.1 then it satisfies
approximately the optimality condition associated to 3.2.

######  Remark 4.2.

There is still another important consequence of the optimality conditions of
3.1. In [45] it was proved that for finite horizon problems the optimality
condition of this method becomes sufficient if the value function is $C^{2}$.
In the following we will see that a similar situation holds for 3.1. Let us
consider $T$ tending to infinity. We formally set $\mathcal{V}_{\infty}$ as
the limit of $\mathcal{V}_{T}$ and $\mathcal{O}_{\infty}$ as the superior
limit of $\mathcal{O}_{T}$ as $T$ tends to infinity.

Let $\theta^{*}\in\mathcal{O}_{\infty}$ be such that it satisfies (4.18) and
set $y(t;y_{0})=y(t;y_{0},u_{v(\theta^{*})})$ for
$(t,y_{0})\in(0,\infty)\times\omega$ . For $T>0$, by (A.7) and (A.5) we deduce
that

$\begin{array}[]{l}\displaystyle\ell(y(t;y_{0}))+\frac{1}{2\beta}|B^{\top}\nabla
v(\theta^{*})(y(t;y_{0}))|^{2}+\nabla\mathcal{V}_{T}(y(t;y_{0}))^{\top}\cdot(f(y(t;y_{0}))-\frac{1}{\beta}BB^{\top}\nabla
v(\theta^{*})(y(t;y_{0})))\\\
\displaystyle=\ell(y(T;y_{0}))+\frac{1}{2\beta}|B^{\top}\nabla
v(\theta^{*})(y(T;y_{0}))|^{2}.\end{array}$

Further assuming that

$\lim_{t\to\infty}\ell(y(t;y_{0}))+\frac{1}{2\beta}|B^{\top}\nabla
v(\theta^{*})(y(t;y_{0}))|^{2}=0,$ (4.20)

and letting $T$ tend to infinity we obtain

$\ell(y(t;y_{0}))+\frac{1}{2\beta}|B^{\top}\nabla
v(\theta^{*})(y(t;y_{0}))|^{2}+\nabla\mathcal{V}_{\infty}(y(t;y_{0}),u_{v(\theta^{*})})^{\top}\cdot(f(y(t;y_{0}))-\frac{1}{\beta}BB^{\top}\nabla
v(\theta^{*})(y(t;y_{0})))=0.$ (4.21)

If there exists $\hat{\vartheta}\in\Theta$ such that
$D_{\theta}v(\theta^{*})\cdot\hat{\vartheta}=v(\theta^{*})-\mathcal{V}_{\infty}(\cdot,u_{v(\theta^{*})})$,
then plugging $\hat{\vartheta}$ in (4.18) leads to

$B^{\top}\nabla
v(\theta^{*})(y(t;y_{0}))=B^{\top}\nabla\mathcal{V}_{\infty}(y(t;y_{0}),u_{v(\theta^{*})})\mbox{
for all }(t,y_{0})\in(0,\infty)\times\omega.$ (4.22)

We note that if $D_{\theta}v(\theta^{*})$ is surjective then the existence of
such a $\hat{\vartheta}$ is equivalent to assuming that
$\mathcal{V}_{\infty}(\cdot,u_{v(\theta^{*})})\in\Theta$. Combining (4.22) and
(4.21) we obtain

$\begin{array}[]{l}\ell(y(t;y_{0}))+\frac{1}{2\beta}|B^{\top}\nabla\mathcal{V}_{\infty}(y(t;y_{0}),u_{v(\theta^{*})})|^{2}\\\\[6.45831pt]
\quad+\nabla\mathcal{V}_{\infty}(y(t;y_{0}),u_{v(\theta^{*})})^{\top}\cdot(f(y(t;y_{0}))-\frac{1}{\beta}BB^{\top}\nabla\mathcal{V}_{\infty}(y(t;y_{0}),u_{v(\theta^{*})}))=0,\end{array}$
(4.23)

which implies that $\mathcal{V}_{\infty}(\cdot,u_{v(\theta^{*})})$ satisfies
equation (2.4) on the set
$\\{y(t;y_{0})\,|\,(t,y_{0})\in(0,\infty)\times\omega\\}$. Further by the
verification theorem (see Theorem I.7.1 in [50]), if $\omega$ is large enough
in the sense that
$\\{y(t;y_{0})\,|\,(t,y_{0})\in(0,\infty)\times\omega\\}\supset\\{y^{*}(t;y_{0})\,|\,(t,y_{0})\in(0,\infty)\times\omega\\}$,
where $y^{*}$ is optimal for (2.1), then
$\mathcal{V}_{\infty}(\cdot,u_{v(\theta^{*})})=V$. In this case
$u_{v(\theta^{*})}$ provides an optimal feedback law for (2.1) and it is an
optimal solution of 3.1.

##  5 Obstacle problem

In this section we introduce a parameterized family of control problems with
different degrees of smoothness of the value function depending on a penalty
parameter. This example was already introduced in [46], for that reason we
will follow the same definition.

For $\gamma\in[0,\infty)$, we consider the control problem

$\min_{\begin{array}[]{c}u\in L^{2}((0,\infty);\mathbb{R}^{2}),\\\
y^{\prime}=u,\
y(0)=y_{0},\end{array}}\int_{0}^{\infty}\ell_{\gamma}(y(t))dt+\frac{\beta}{2}\int_{0}^{\infty}|u(t)|^{2}dt,$
(5.1)

where the running cost for the state variable is given by

$\ell_{\gamma}(y)=\frac{1}{2}|y|^{2}\left(1+\gamma\psi\left(\frac{|y-z|}{\sigma}\right)\right),\text{
with
}\psi(s)=\left\\{\begin{array}[]{ll}\exp\left(-\frac{1}{1-s^{2}}\right)&\mbox{
if }|s|<1\\\ 0&\mbox{ if }|s|\geqslant 1,\end{array}\right.$

and $z\in\mathbb{R}^{2}$ satisfies $z_{1}<-\sigma$, $z_{2}=0$, for some
$\sigma>0$.

The value function of this problem is denoted by $V_{\gamma}$. In [46], local
Lipschitz regularity of $V_{\gamma}$ for all $\gamma>0$ and its non-
differentiability for $\gamma$ sufficiently large were proved. Here we will
prove that for $\gamma$ close to $0$ the function $V_{\gamma}$ is
$C^{2}(\mathbb{R}^{2})$. That is, we will establish that for
$\gamma\in[0,\gamma_{s})$ the function $V_{\gamma}$ is
$C^{2}(\mathbb{R}^{2})$, where

$\gamma_{s}=\frac{1}{\sup_{y\in
B(z,\sigma)}\left\\{2|y|\cdot|\nabla\tilde{\psi}(y)|+\frac{1}{2}|y|\cdot|\nabla^{2}\tilde{\psi}(y)|\right\\}}$
(5.2)

with $\tilde{\psi}(y)=\psi\left(\frac{|y-z|}{\sigma}\right)$.

######  Proposition 5.1.

Let $\gamma\in[0,\gamma_{s})$. The value function $V_{\gamma}$ of (5.1) is
$C^{\infty}(\mathbb{R}^{2})$.

###### Proof.

By [46] we already know that for all $y_{0}\in\mathbb{R}^{d}$ and
$\gamma\in[0,\infty)$ problem (5.1) has at least one solution and $V_{\gamma}$
is locally Lipschitz continuous. Additionally, we observe that

$\nabla^{2}\ell_{\gamma}(y)=I(1+\gamma\tilde{\psi}(y))+2\gamma
y\otimes\nabla\tilde{\psi}(y)+\gamma\frac{|y|^{2}}{2}\nabla^{2}\tilde{\psi}(y)$

where $I$ is the $2\times 2$ identity matrix. Multiplying both sides of the
previous equality by $x\in\mathbb{R}^{2}\setminus{(0,0)}$ we see that

$x^{\top}\cdot\nabla^{2}\ell_{\gamma}(y)\cdot
x\geqslant\left((1+\gamma\tilde{\psi}(y))-2\gamma|y|\cdot|\nabla\tilde{\psi}(y)|-\gamma\frac{|y|^{2}}{2}\cdot|\nabla^{2}\tilde{\psi}(y)|\right)|x|^{2}.$

Using that the support of $\tilde{\psi}$ is $B(z,\sigma)$ in the previous
inequality, we obtain for $\gamma\in[0,\gamma_{s})$ that

$x^{\top}\cdot\nabla^{2}\ell_{\gamma}(y)\cdot
x\geqslant\left(1-\gamma\sup_{y\in
B(z,\sigma)}\left\\{2|y|\cdot|\nabla\tilde{\psi}(y)|+\frac{1}{2}|y|\cdot|\nabla^{2}\tilde{\psi}(y)|\right\\}\right)|x|^{2}>0,$
(5.3)

which implies that $\ell_{\gamma}$ is a strictly convex function. Due to the
linearity of the state equation of (5.1), we also get that the objective
function of (5.1) is convex and consequently the optimal control of (5.1) is
unique for $\gamma\in[0,\gamma_{s})$.

In the following, we will use (5.3) to prove that $V_{\gamma}$ is
$C^{2}(\mathbb{R}^{2})$ by means of the implicit function theorem. To this
end, we define set $X=H^{1}((0,\infty);\mathbb{R}^{2})\times
L^{2}((0,\infty);\mathbb{R}^{2})\times
H^{1}((0,\infty);\mathbb{R}^{2})\times$,
$Y=L^{2}((0,\infty);\mathbb{R}^{2})\times
L^{2}((0,\infty);\mathbb{R}^{2})\times\mathbb{R}^{2}\times
L^{2}((0,\infty);\mathbb{R}^{2})$ and a function
$\Phi:X\times\mathbb{R}^{2}\mapsto Y$ by

$\Phi(y,u,p,y_{0})=\left(\begin{array}[]{l}-p^{\prime}+\nabla\ell_{\gamma}(y)\\\
y^{\prime}-u\\\ y(0)-y_{0}\\\ \beta u-p\end{array}\right)$ (5.4)

for $((y,u,p),y_{0})\in X$. We note that $\Phi(y,u,p,y_{0})=0$ is the
optimality condition for (5.1). Clearly the function $\Phi$ is of class
$C^{\infty}$ and its differential with respect to its first three variables is
given by

$D_{y,u,p}(y,u,p,y_{0})\cdot(\delta_{y},\delta_{u},\delta_{p})=\left(\begin{array}[]{l}-\delta_{p}^{\prime}+\nabla^{2}\ell_{\gamma}(y)\delta_{y}\\\
\delta_{y}^{\prime}-\delta_{u}\\\ \delta_{y}(0)\\\
\beta\delta_{u}-\delta_{p}\end{array}\right)$

for $((y,u,p),y_{0})\in X\times\mathbb{R}^{2}$ and
$(\delta_{y},\delta_{u},\delta_{p})\in X$. In order to use the implicit
function theorem we need to prove that for $u^{*}$ the optimal solution of
(5.1), $y^{*}$ its optimal trajectory and $p^{*}$ the respective adjoint
state, the partial sub-differential $D_{y,u,p}(y^{*},u^{*},p^{*})$ is
invertible , i.e., for all $(z_{p},z_{y},z_{y_{0}},z_{u})\in Y$ there exist
only one $(\delta_{y},\delta_{u},\delta_{p})\in X$ such that

$D_{y,u,p}(y^{*},u^{*},p^{*})\cdot(\delta_{y},\delta_{u},\delta_{p})=(z_{p},z_{y},z_{y_{0}},z_{u}).$

We note that this equation is equivalent to the optimality conditions of the
following linear-quadratic problem

$\min_{\begin{array}[]{c}\delta_{u}\in L^{2}((0,\infty);\mathbb{R}^{2}),\\\
\delta_{y}^{\prime}=\delta_{u}+z_{y},\
\delta_{y}(0)=z_{y_{0}}\end{array}}\int_{0}^{\infty}\left(\frac{1}{2}\delta_{y}^{\top}(t)\nabla^{2}\ell_{\gamma}(y^{*}(t))\delta_{y}(t)-z_{p}(t)\cdot\delta_{y}(t)+\frac{\beta}{2}|\delta_{u}(t)|^{2}-\delta_{u}(t)\cdot
z_{u}(t)\right)dt.$ (5.5)

Consider $\tilde{u}$ the optimal solution of (5.1) with initial condition
$z_{y_{0}}$ and $\tilde{y}$ the corresponding trajectory. We note that using
$\delta_{u}=\tilde{u}-z_{y}$ we obtain that $\delta_{y}=\tilde{y}$ satisfies
$\delta_{y}^{\prime}-\delta_{u}=z_{y}$ and

$\left|\int_{0}^{\infty}\left(\frac{1}{2}\delta_{y}^{\top}(t)\nabla^{2}\ell_{\gamma}(y(t))\delta_{y}(t)-z_{p}(t)\cdot\delta_{y}(t)+\frac{\beta}{2}|\delta_{u}(t)|^{2}-\delta_{u}(t)\cdot
z_{u}(t)\right)dt\right|<+\infty.$

This together with the positive definiteness of $\nabla^{2}\ell_{\gamma}$ and
the weak lower semi-continuity of the objective function imply that (5.5) has
a unique optimal solution. Then by the implicit function theorem there exists
an open neighborhood $\Gamma\subset\mathbb{R}^{2}$ of $y_{0}$ and functions
$\hat{u}:\Gamma\mapsto L^{2}((0,\infty);\mathbb{R}^{2})$,
$\hat{y}:\Gamma\mapsto H^{1}((0,\infty);\mathbb{R}^{2})$ and
$\hat{p}:\Gamma\mapsto H^{1}((0,\infty);\mathbb{R}^{2})$ of class
$C^{\infty}$, such that for all $z_{0}\in\Gamma$ we have that
$\Phi(y(z_{0}),u(z_{0}),p(z_{0}),z_{0})=0$.Since for $\gamma\in[0,\gamma_{s})$
problem (5.1) is strictly convex, we have that $u(z_{0})$ is its optimal
solution, $y(z_{0})$ its optimal trajectory, and $p(z_{0})$ the corresponding
adjoint state. Hence,
$V_{\gamma}(z_{0})=\int_{0}^{\infty}\left(\ell_{\gamma}(y(z_{0}))+\frac{\beta}{2}|u(z_{0})|^{2}\right)dt$
and by the chain rule we obtain that $V_{\gamma}$ is $C^{\infty}$ in $\Gamma$,
Since $y_{0}$ is arbitrary we deduce that $V_{\gamma}$ is
$C^{\infty}(\mathbb{R}^{2})$. ∎

As mentioned before, in [46] it was proved that for $\gamma$ large enough the
value function $V_{\gamma}$ is non differentiable. Specifically, for each such
$\gamma$ there exists $x_{\gamma}$ such that $V_{\gamma}$ is non-
differentiable in $\\{(x,0)\in\mathbb{R}^{2}:x>x_{\gamma}\\}$. However, it was
also proved in [46] that for all $\gamma>0$ 3.1 and 3.2 hold. Consequently, by
using a sequence of setting satisfying 3.3 we obtain by means of 4.1 that
convergence in the sense of (3.7) holds for 3.1. On the other hand, for
Methods 3.2 and 3.3 and $\gamma\in[0,\gamma_{s})$ we obtain that 3.4 holds
provided the sequence of settings is chosen appropriately.

This together with 5.1 allows us to utilize problem (5.1) to numerically
investigate the behavior of Methods 3.1 \- 3.3 for different choices of
$\gamma$, which account for different degrees of smoothness of $V_{\gamma}$.

We close this section by establishing another important property of the value
function.

######  Proposition 5.2.

The value function associated to (5.1) is semiconcave on every compact subset
of $\mathbb{R}^{2}$.

###### Proof.

Let $C\subset\mathbb{R}^{2}$ be compact and denote by
$(y^{*}(\cdot;y_{0},u^{*}(\cdot;y_{0})))$ the optimal state-control pairs with
$y_{0}\in C$. Since the optimal trajectories associated to (5.1) are globally
exponentially stable there exists $\hat{T}$ such that
$y^{*}(t;y_{0})<|z|-\sigma$ for all $(t,y_{0})\in[\hat{T},\infty)\times C$,
and consequently
$V_{\gamma}(y_{0})=V_{0}(y_{0})=\frac{\sqrt{\beta}}{2}|y_{0}|^{2}$, [46, Lemma
9.2, Lemma 9.3]. Here $B(0,|z|-\sigma)$ denotes the open ball with center at
the origin and radius $|z|-\sigma$. Further
$|u^{*}(t;y_{0})|=\frac{2}{\beta}\ell_{\gamma}(y^{*}(t;y_{0}))$, [46,
Proposition 9.1], and thus for a constant $K$ we have
$|u^{*}(t;y_{0})|\leqslant K$, for all $(t,y_{0})\in[0,\infty)\times C$.

Let us consider

$\min_{\begin{array}[]{c}u\in
L^{2}((0,\hat{T});\mathbb{R}^{2}),\,|u(t)|\leqslant K+1\\\ y^{\prime}=u,\
y(0)=y_{0}\in
C\end{array}}\int_{0}^{\hat{T}}\ell_{\gamma}(y(t))dt+\frac{\beta}{2}\int_{0}^{\hat{T}}|u(t)|^{2}dt+\hat{V}_{0}(y(\hat{T},y_{0})),$
(5.6)

where $\hat{V}_{0}$ denotes the extension of $V_{0}$ outside of
$B(0,|z|-\sigma)$ by $\hat{V}_{0}(y)=\frac{\sqrt{\beta}}{2}|y|^{2}$.

For $y_{0}\in C$ this problem is equivalent to (5.1). Here we use Bellman’s
principle, the fact that $y^{*}(\hat{T};y_{0})\in B(0,|z|-\sigma)$, and that
$|u^{*}|_{L^{\infty}(0,\infty)}\leqslant K$, so that the control constraint in
(5.6) is not active at the solutions. The value function restricted to $C$
coincides for these two problems. We next utilize that (5.6) can be considered
as finite horizon problem with terminal penalty given by $\hat{V}_{0}$. It
then follows from [51, Theorem 7.4.11] that the value function for (5.6) is
semiconcave on $C$, provided that the assumptions of that theorem are
satisfied. In this respect it suffices to recall that any function which is
$C^{1}(A)$ in an open set $A$ containing $C$ is semiconcave on $C$, [51,
Proposition 2.1.2]. This is directly applicable for $\ell\in
C^{1}(\mathbb{R}^{2})$ as well as $\hat{V}_{0}\in C^{\infty}(\mathbb{R}^{2})$.
The remaining assumptions on the dynamical system and the running cost are
trivially satisfied. ∎

##  6 Numerical experiments

Numerical tests for Methods 3.1 and 3.2 are available in the literature, we
refer to eg. [16, 17, 22, 46, 48, 14, 47]. Thus the strength and validity of
them for practical realisation has already been established. Here we make an
effort to compare these methods with respect to the influence of the
regularity of the value function. For this purpose, we employ Equation 5.1 and
recall that the parameter $\gamma$ acts as a control for the regularity of
$V_{\gamma}$. It would be of interest to also compare these methods for high
dimensional problems. This, however, is not within the scope of this research.

We consider problem 5.1 with $z_{1}=-2$, $\sigma=1$, and set
$\omega=(-5,5)\times(-2,2)$. Since the dimensionality is low, a regular grid
is used to approximate the integrals appearing in Methods 3.1 \- 3.3. For high
dimensional problems this is not feasible and other discretization method
should be used instead, as for example Monte-Carlo integration, which is the
usual choice for machine learning methods and has been used in [48, 45]. We
report on the behavior of the obtained feedback laws for each of the methods
with $\gamma\in\\{10^{-3},10^{-2},10^{-1},1,10,10^{2},10^{3}\\}$.

In the following we further describe the setting of the methods. For the
Banach spaces $\Theta$ and the parametrization we follow 3.2. Thus the setting
depends on two parameters, which are the order $k$ and the degree $n$. Here we
consider $k\in\\{1,2\\}$ and $n\in\\{10,20,30,40\\}$. Recall that the choice
of $k$ connects the smoothness of the value function. The penalty function of
the setting is
$\mathcal{P}(\theta)=\frac{1}{2}\left(\frac{1}{2}|\theta|_{2}^{2}+|\theta|_{1}\right)$
and the our tests were carried out with penalty parameter $\alpha$ chosen from
the set $\\{10^{-1},10^{-3},10^{-5},10^{-7},10^{-9}\\}$. In all the cases we
consider the time horizon with $T=1$. Due to the amount of data we only
present for each problem, order $k$, and $\gamma$ the degree $n$ and the
penalty parameter of the setting $\alpha$ with the best $NMAE_{V}$ in the
training set. We present this information in Table 1, where it can be seen
that typically the best performance is achieved for $\alpha$ around $10^{-}7$
and with the degree equal to 30 or 40.

$\gamma$ | Problem | Order | Degree | $\alpha$ |  $NMAE_{V}$ training[%] |  $NMRSE_{c}$ training[%] |  Training Time [s]
---|---|---|---|---|---|---|---
1.0e-03 | Problem 3.1 | 1 | 10 | 1.0e-09 | 3.48e-04 | 3.32e-01 | 3.64e+02
1.0e-03 | Problem 3.1 | 2 | 30 | 1.0e-05 | 4.93e-03 | 7.13e-01 | 5.47e+02
1.0e-03 | Problem 3.2 | 1 | 10 | 1.0e-09 | 2.43e-04 | 1.90e-01 | 2.20e+02
1.0e-03 | Problem 3.2 | 2 | 10 | 1.0e-07 | 1.85e-03 | 4.20e-01 | 2.08e+02
1.0e-03 | Problem 3.3 | 1 | 10 | 1.0e-03 | 1.67e-03 | 3.77e-01 | 2.30e+00
1.0e-03 | Problem 3.3 | 2 | 40 | 1.0e-03 | 1.06e-03 | 3.02e-01 | 5.37e+00
1.0e-02 | Problem 3.1 | 1 | 10 | 1.0e-09 | 3.31e-04 | 3.07e-01 | 3.88e+02
1.0e-02 | Problem 3.1 | 2 | 30 | 1.0e-05 | 4.90e-03 | 7.01e-01 | 5.78e+02
1.0e-02 | Problem 3.2 | 1 | 10 | 1.0e-09 | 2.26e-04 | 1.79e-01 | 2.29e+02
1.0e-02 | Problem 3.2 | 2 | 10 | 1.0e-07 | 1.85e-03 | 4.11e-01 | 2.18e+02
1.0e-02 | Problem 3.3 | 1 | 10 | 1.0e-03 | 1.50e-03 | 3.53e-01 | 3.05e+00
1.0e-02 | Problem 3.3 | 2 | 40 | 1.0e-03 | 1.00e-03 | 2.88e-01 | 6.47e+00
1.0e-01 | Problem 3.1 | 1 | 10 | 1.0e-09 | 1.06e-03 | 3.82e-01 | 4.37e+02
1.0e-01 | Problem 3.1 | 2 | 30 | 1.0e-05 | 4.62e-03 | 6.72e-01 | 7.10e+02
1.0e-01 | Problem 3.2 | 1 | 10 | 1.0e-05 | 5.79e-04 | 2.61e-01 | 2.72e+02
1.0e-01 | Problem 3.2 | 2 | 10 | 1.0e-07 | 2.68e-03 | 5.13e-01 | 2.39e+02
1.0e-01 | Problem 3.3 | 1 | 10 | 1.0e-03 | 2.06e-03 | 4.19e-01 | 4.54e+00
1.0e-01 | Problem 3.3 | 2 | 20 | 1.0e-09 | 2.89e-03 | 4.97e-01 | 6.53e+00
1.0e+00 | Problem 3.1 | 1 | 40 | 1.0e-07 | 4.41e-03 | 8.20e-01 | 1.70e+03
1.0e+00 | Problem 3.1 | 2 | 20 | 1.0e-05 | 3.88e-02 | 2.07e+00 | 9.69e+02
1.0e+00 | Problem 3.2 | 1 | 40 | 1.0e-07 | 2.98e-03 | 6.39e-01 | 1.68e+03
1.0e+00 | Problem 3.2 | 2 | 40 | 1.0e-05 | 4.56e-02 | 2.15e+00 | 7.09e+02
1.0e+00 | Problem 3.3 | 1 | 30 | 1.0e-03 | 1.17e-02 | 9.40e-01 | 8.71e+01
1.0e+00 | Problem 3.3 | 2 | 40 | 1.0e-07 | 1.72e-02 | 1.25e+00 | 3.11e+01
1.0e+01 | Problem 3.1 | 1 | 40 | 1.0e-09 | 1.05e-01 | 3.04e+00 | 3.49e+03
1.0e+01 | Problem 3.1 | 2 | 30 | 1.0e-09 | 8.23e-01 | 6.48e+00 | 3.42e+03
1.0e+01 | Problem 3.2 | 1 | 30 | 1.0e-07 | 1.55e-01 | 2.97e+00 | 3.10e+03
1.0e+01 | Problem 3.2 | 2 | 40 | 1.0e-09 | 1.65e+00 | 5.98e+00 | 3.04e+03
1.0e+01 | Problem 3.3 | 1 | 30 | 1.0e-03 | 6.72e-01 | 9.01e+00 | 1.28e+02
1.0e+01 | Problem 3.3 | 2 | 30 | 1.0e-07 | 4.25e-01 | 6.04e+00 | 1.93e+02
1.0e+02 | Problem 3.1 | 1 | 40 | 1.0e-07 | 7.09e-01 | 6.49e+00 | 5.66e+03
1.0e+02 | Problem 3.1 | 2 | 20 | 1.0e-09 | 2.98e+00 | 1.14e+01 | 5.36e+03
1.0e+02 | Problem 3.2 | 1 | 40 | 1.0e-05 | 1.18e+00 | 5.07e+00 | 6.01e+03
1.0e+02 | Problem 3.2 | 2 | 40 | 1.0e-09 | 2.02e+01 | 7.90e+00 | 3.93e+03
1.0e+02 | Problem 3.3 | 1 | 10 | 1.0e-01 | 8.74e+00 | 2.48e+01 | 7.06e+00
1.0e+02 | Problem 3.3 | 2 | 10 | 1.0e-09 | 9.64e+00 | 3.14e+01 | 9.14e+01
1.0e+03 | Problem 3.1 | 1 | 40 | 1.0e-09 | 2.66e+00 | 9.34e+00 | 8.16e+03
1.0e+03 | Problem 3.1 | 2 | 40 | 1.0e-07 | 8.31e+00 | 2.05e+01 | 2.19e+04
1.0e+03 | Problem 3.2 | 1 | 40 | 1.0e-05 | 9.05e+00 | 7.38e+00 | 9.96e+03
1.0e+03 | Problem 3.2 | 2 | 20 | 1.0e-07 | 1.93e+02 | 9.25e+00 | 2.35e+03
1.0e+03 | Problem 3.3 | 1 | 10 | 1.0e-01 | 2.28e+01 | 4.31e+01 | 9.15e+00
1.0e+03 | Problem 3.3 | 2 | 10 | 1.0e-09 | 2.43e+01 | 4.56e+01 | 1.12e+02
Table 1: Each row represents a the degree and $\alpha$ with the lowest
training NMAE for the correspoing problem and order.

The integrals over $\omega$ appearing in the learning problems are
approximated over a regular grid of size $16\times 16$. We call these points
the training set. In order two study the performance of the feedback produced
by each method we also consider a finer regular grid of size $32\times 32$,
which is the test set. In both the training and test sets the performance of
each method is evaluated by the following error measures:

$NMAE_{V}=\sum_{i=1}^{N}|V(y_{0}^{i})-\mathcal{V}_{T}(y_{0}^{i},v(\theta))|\cdot\frac{1}{\sum_{i=1}^{N}V(y_{0}^{i})},$
(6.1) $NMRSE_{c}=\sum_{i=1}^{N}\left\lVert
u_{i}^{*}-\hat{u}_{i}\right\rVert_{L^{2}((0,T);\mathbb{R}^{2})}\cdot\frac{1}{\sum_{i=1}^{N}\left\lVert
u_{i}^{*}\right\rVert_{L^{2}((0,T);\mathbb{R}^{2})}},$ (6.2)

where $\\{y_{0}^{i}\\}_{i=1}^{n}$ is a set of initial conditions which in this
case could be either the training or the test sets,
$\\{u_{i}^{*}\\}_{i=1}^{N}$ are the optimal controls and
$\\{\hat{u}_{i}\\}_{i=1}^{N}$ are the controls obtained by applying the
respective feedback. For obtaining a reference value for the optimal value
function $V$ and the optimal controls for the chosen initial conditions we
numerically solve the discretized open loop control problem. For this purpose
the $ODEs$ and the integral of the running cost are discretized by an explicit
Euler scheme with time step $\frac{1}{400}$. In accordance with Theorems 4.1
and 4.2 we expect that for small values of $\gamma$ all three approaches will
achieve a good performance since in this case $V_{\gamma}$ is smooth. On the
other hand, $V_{\gamma}$ is non-smooth for large values of $\gamma$, and only
Theorem 4.1 is applicable, which implies that we only have a guarantee of good
performance for Method 3.1.

For solving the learning problems we employ a proximal gradient method
together with the Barzilai-Borwein step size and a non-linear backtracking
line search. This algorithm is described in [48] and [45] in the context of
3.1. Following the notation of Algorithm 1 in [45] we use $k=3$,
$tol=10^{-6}$, $\kappa=10^{-3}$ and $\xi=\frac{1}{2}$, with the same
tolerances measures as in [45]. Concerning initialization, for each setting we
first solve with $\gamma=10^{-3}$ using as initial guess $\theta^{0}=0$, and
subsequently for $\gamma>10^{-3}$ the solution corresponding to the previous
$\gamma$ value is used as initial guess.

(a)

(b)

(c)

Figure 1: Level sets of $V_{\gamma}$, the color of the level sets are in
logarithmic scale.

We start the discussion of the experiments with the level sets of the value
function for selected values of regularisation parameters
$\gamma=10^{-3},0,10^{3}$. This is relevant to understand the distribution of
the errors for each method which will be reported later on. In Figure 1 the
level sets for selected $\gamma$-values together with the boundary of the
obstacle are presented. We observe that on the right hand side of the $y-$axis
the level sets do not change with $\gamma$. This is concordant with the fact
that in this region the obstacle has no effect. On the other hand, on the left
hand side we observe that in a neighborhood of the $x-$axis the level sets
vary from being convex to non-convex. This is consistent with the fact that on
this axis the value function is non-differentiable for large $\gamma$.
Further, for $\gamma=10^{3}$ the value function has a local maximum close to
the center of the obstacle. Indeed, since $V_{\gamma}$ converges to
$V_{\infty}$ it is expected that $V_{\gamma}$ blows up inside the obstacle as
$\gamma$ tends to infinity.

(a)

(b)

(c)

(d)

(e)

Figure 2: Errors and training times. In all the figures both axes are in
logarithmic scale.

(a)

(b)

(c)

(d)

Figure 3: Errors normalized by the error of Method 3.1. In all the figures
both axes are in logarithmic scale.

To compare the performance graphically, in Figure 2 the errors and the
training time for each method are depicted. We observe that typically $k=1$
achieves a better performance than $k=2$, especially for $\gamma\geqslant 10$.
Additionally, we note that the error increases with $\gamma$. Methods 3.1 and
3.2 with $k=1$ achieved the best performance in both error measures. To take a
closer look we turn our attention to Figure 3 where the errors are normalized
by dividing them by the error obtained by using Method 3.1 with k=1. From this
we see that for $\gamma\leqslant 1$ the method with the best performance in
terms of $NMAE_{V}$ and $NMRSE_{c}$ control is Method 3.2 with $k=1$. In
contrast, for $\gamma>1$, Method 3.1 with $k=1$ has the best $NMAE_{V}$ and
moreover the distance with respect to Method 3.2 increases. We point out that
also for $\gamma>1$ the fact that Method 3.2 has best $NMRSE_{c}$ control does
not imply the same with respect to $NMAE_{V}$. Concerning the efficiency of
the methods in Figure 2(2(e)) we provide the trainings time for each
experiment, that is, the time that it takes to achieve the stooping criterion.
In this figure we see that Method 3.3 is the one with the shortest training
time. For Methods 3.1 and 3.2 the training times tend to increase with
$\gamma$. This can be a consequence of decreasing regularity and the fact that
the value function blows up inside the obstacle as $\gamma$ tends to infinity.

(a)

(b)

(c)

(d)

(e)

(f)

Figure 4: Scatter plots for $\gamma=10^{3}$.

(a)

(b)

(c)

(d)

(e)

(f)

Figure 5: Error distribution for $\gamma=10^{-3}$.

(a)

(b)

(c)

(d)

(e)

(f)

Figure 6: Error distribution for $\gamma=1$.

(a)

(b)

(c)

(d)

(e)

(f)

Figure 7: Error distribution for $\gamma=10^{3}$.

The previous analysis does not address the spatial distribution of the error.
This information is given in Figures 5-7 where the distribution of the error
between the value functions $V_{\gamma}$ and
$\mathcal{V}_{T}(\cdot,\theta^{*})$ is depicted, with $\theta^{*}$ the
solution obtained by the respective experiment. From Figure 5 we observe that
the error is uniformly smaller than $10^{-3}$ for $\gamma\leqslant 10^{-3}$,
and we can report that this is also the case as long as $\gamma\leqslant 1$.
Starting from $\gamma=10^{1}$ we see from Figures 6-7 that the error in the
region left to the $y-$axis tends to increase with $\gamma$. Moreover, the
error is higher for the points close to the $x-$axis which are to the left of
the $y-$axis and also inside the obstacle than those outside. As explained
before, this can be attributed to the fact that the value function is non
smooth on the $x-axis$ for all $x$ left of the obstacle and that it blows up
inside the obstacle as $\gamma$ tends to infinity.

Figures 4 show the scatter plot and the dispersion between $V_{\gamma}$ and
$\mathcal{V}_{T}(\cdot,\theta^{*})$ evaluated at the test points for
$\gamma=10^{3}$. In all the cases the points in the interior of the obstacle
tend to be further from the identity line than those outside of the obstacle.
For Method 3.1 this effect is weaker than for the others.

##  7 Conclusions

In this work convergence in the sense of (3.7) is analyzed for three machine
learning approaches. They differ with respect to the assumption on the
regularity of the value function. Further, now that more techniques to compute
suboptimal feedback controls guided by the Bellman principle and hence the HJB
equation, become available, it is of interest to start towards the challenging
goal of comparing methods. For this purpose we present a family of control
problems in Section 5 indexed by a penalty parameter $\gamma>0$ for which the
value function is smooth if $\gamma$ close to 0, and while it is non-
differentiable but semi-concave for $\gamma$ large. This family of control
problem was used to compare the performance of the three methods. From these
experiments it was observed that for $\gamma$ close to 0, which correspond to
the case of smooth value functions, Method 3.2 has a better performance than
the other approaches, but for larger $\gamma$ Method 3.1 is the one with a
better performance. This is coherent with the convergence results.

Further extensions of the work presented here are possible. For instance, the
discretized version of the learning problems was not studied. This is
particularly relevant for high-dimensional problems, since in this case the
discretization is usually carried out by using Monte-Carlo integration.
Problems with constraints on either the state or the control can be considered
by first using an exterior penalty method as was used for the construction of
problem (5.1). Finally, the extension of this problem to higher dimension is
of interest, because it could help to study the behavior of the learning
problems in the context of high-dimensions and non-smoothness of the value
function.

## References

  * [1] F. Bonnans, P. Chartier, H. Zidani, Discrete approximation for a class of the Hamilton-Jacobi equation for an optimal control problem of a differential-algebraic system, Control and Cybernetics 32 (1) (2003) 33–55.
  * [2] M. Falcone, R. Ferretti, Semi-Lagrangian approximation schemes for linear and Hamilton-Jacobi equations, SIAM, Philadelphia, 2013.
  * [3] A. Alla, M. Falcone, D. Kalise, An efficient policy iteration algorithm for dynamic programming equations, SIAM J. Control Optim. 37 (1) (2015) A181–A200.
  * [4] R. W. Beard, G. N. Saridis, J. T. Wen, Galerkin approximation of the generalized Hamilton-Jacobi-Bellman equation, Automatica 33 (12) (1997) 2159–2177.
  * [5] M. Puterman, S. Brumelle, On the convergence of policy iteration in stationary dynamic programming, Math. Oper. Res. 4 (1) (1979) 60–69.
  * [6] M. Santos, J. Rust, Convergence properties of policy iteration, SIAM J. Control Optim. 42 (6) (2004) 2094–2115.
  * [7] Y. T. Chow, J. Darbon, S. Osher, W. Yin, Algorithm for overcoming the curse of dimensionality for time-dependent non-convex Hamilton– Jacobi equations arising from optimal control and differential games problems, J. Sci. Comput. 73 (2017) 617–643.
  * [8] Y. T. Chow, J. Darbon, S. Osher, W. Yin, Algorithm for overcoming the curse of dimensionality for state-dependent Hamilton–Jacobi equations, J. Comput. Phys 387 (2019) 376–409.
  * [9] Y. T. Chow, W. Li, S. Osher, W. Yin, Algorithm for Hamilton–Jacobi equations in density space via a generalized Hopf formula, J. Sci. Comput. 80 (2019) 1195–1239.
  * [10] J. Darbon, S. Osher, Algorithms for overcoming the curse of dimensionality for certain Hamilton–Jacobi equations arising in control theory and elsewhere, Res. Math. Sci. 3 (2016) 19.
  * [11] J. Han, A. Jentzen, E. Weinan, Solving high-dimensional partial differential equations using deep learning, Proc. Nat. Acad. Sci. USA 115 (34) (2018) 8505–8510.
  * [12] J. Darbon, G. P. Langlois, T. Meng, Overcoming the curse of dimensionality for some Hamilton–Jacobi partial differential equations via neural network architectures, Res. Math. Sci. 7 (2020) 1–50.
  * [13] N. Nüsken, L. Richter, Solving high-dimensional Hamilton–Jacobi– Bellman pdes using neural networks: Perspectives from the theory of controlled diffusions and measures on path space (2020).
URL https://arxiv.org/abs/2005.05409

  * [14] D. Onken, L. Nurbekyan, X. Li, S. W. Fung, S. Osher, L. Ruthotto, A neural network approach applied to multi-agent optimal control, in: 2021 European Control Conference (ECC), 2021, pp. 1036–1041.
  * [15] K. Ito, C. Reisinger, Y. Zhang, A neural network-based policy iteration algorithm with global h2-superlinear convergence for stochastic games on domains, Found. Comput. Math. 21 (2) (2021) 331–374.
  * [16] K. Kunisch, D. Walter, Semiglobal optimal feedback stabilization of autonomous systems via deep neural network approximation, ESAIM - Control Optim. Calc. Var. 27 (2021) 16.
  * [17] L. Ruthotto, S. J. Osher, W. Li, L. Nurbekyan, S. W. Fung, A machine learning framework for solving high-dimensional mean field game and mean field control problems, Proc. Natl. Acad. Sci. 117 (17) (2020) 9183–9193.
  * [18] G. Chen, Deep neural network approximations for the stable manifolds of the hamilton-jacobi equations (2023).  arXiv:2007.15350.
  * [19] Y. Zhao, J. Han, Offline supervised learning v.s. online direct policy optimization: A comparative study and a unified training paradigm for neural network-based optimal feedback control (2023).  arXiv:2211.15930.
  * [20] T. Nakamura-Zimmerer, Q. Gong, W. Kang, Qrnet: Optimal regulator design with LQR-augmented neural networks, EEE Control Syst. Lett. 5 (4) (2021) 1303–1308.
  * [21] T. Nakamura-Zimmerer, Q. Gong, W. Kang, Adaptive deep learning for high-dimensional Hamilton–Jacobi–Bellman equations, SIAM J. Sci. Comput. 43 (2) (2021) A1221–A1247.
  * [22] B. Azmi, D. Kalise, K. Kunisch, Optimal feedback law recovery by gradient-augmented sparse polynomial regression, J. Mach. Learn. Res. 22 (48) (2021) 1–32.
  * [23] W. Kang, Q. Gong, T. Nakamura-Zimmerer, Algorithms of data development for deep learning and feedback design, Phys. D: Nonlinear Phenom. 425 (2021) 132955.
  * [24] G. Albi, S. Bicego, D. Kalise, Gradient-augmented supervised learning of optimal feedback laws using state-dependent Riccati equations, IEEE Control Systems Letters 6 (2022) 836–841.
  * [25] T. Ehring, B. Haasdonk, Hermite kernel surrogates for the value function of high-dimensional nonlinear optimal control problems (2023).  arXiv:2305.06122.
  * [26] M. Akian, S. Gaubert, A. Lakhoua, The max-plus finite element method for solving deterministic optimal control problems: Basic properties and convergence analysis, SIAM J. Control Optim. 47 (2) (2008) 817–848.
  * [27] S. Gaubert, W. McEneaney, Z. Qu, Curse of dimensionality reduction in max-plus based approximation methods: Theoretical estimates and improved pruning algorithms, in: 2011 50th IEEE Conference on Decision and Control and European Control Conference, 2011, pp. 1054–1061.
  * [28] P. M. Dower, W. M. McEneaney, H. Zhang, Max-plus fundamental solution semigroups for optimal control problems, in: 2015 Proceedings of the Conference on Control and its Applications, 2015, pp. 368–375.
  * [29] D. Kalise, K. Kunisch, Polynomial approximation of high-dimensional Hamilton–Jacobi–Bellman equations and applications to feedback control of semilinear parabolic pdes, SIAM J. Sci. Comput. 40 (2) (2018) A629–A652.
  * [30] D. Kalise, S. Kundu, K. Kunisch, Robust feedback control of nonlinear pdes by numerical approximation of high-dimensional Hamilton–Jacobi–Isaacs equations, SIAM J. Appl. Dyn. Syst. 19 (2) (2020) 1496–1524.
  * [31] M. B. Horowitz, A. Damle, J. W. Burdick, Linear Hamilton Jacobi Bellman equations in high dimensions, in: 2014 53rd IEEE Conf. Decis. Control, 2014, pp. 5880–5887.
  * [32] E. Stefansson, Y. P. Leong, Sequential alternating least squares for solving high dimensional linear Hamilton–Jacobi–Bellman equation, in: 2016 IEEE/RSJ Int. Conf. Intell. Robots Syst. (IROS), 2016, pp. 3757–3764.
  * [33] A. Gorodetsky, S. Karaman, Y. Marzouk, High-dimensional stochastic optimal control using continuous tensor decompositions, Int. J. Robot. Res. 37 (2-3) (2018) 340–377.
  * [34] S. Dolgov, D. Kalise, K. Kunisch, Tensor decomposition methods for high-dimensional Hamilton–Jacobi–Bellman equations, SIAM J. Sci. Comput. 43 (3) (2021) A1625–A1650.
  * [35] M. Oster, L. Sallandt, R. Schneider, Approximating the stationary Hamilton–Jacobi–Bellman equation by hierarchical tensor products, arXiv (2019).
URL https://arxiv.org/abs/1911.00279

  * [36] M. Oster, L. Sallandt, R. Schneider, Approximating optimal feedback controllers of finite horizon control problems using hierarchical tensor formats, SIAM J. Sci. Comput. 44 (3) (2022) B746–B770.
  * [37] A. Alla, M. Falcone, S. Volkwein, Error analysis for pod approximations of infinite horizon problems via the dynamic programming approach, SIAM Journal on Control and Optimization 55 (5) (2017) 3091–3115.
  * [38] K. Kunisch, S. Volkwein, L. Xie, Hjb-pod-based feedback design for the optimal control of evolution problems, SIAM J. Appl. Dyn. Syst. 3 (4) (2004) 701–722.
  * [39] A. Alla, M. Falcone, L. Saluzzi, An efficient dp algorithm on a tree-structure for finite horizon optimal control problems, SIAM J. Sci. Comput. 41 (4) (2019) A2384–A2406.
  * [40] Bokanowski, J. Garcke, M. Griebel, I. Klompmaker, An adaptive sparse grid semi- Lagrangian scheme for first order Hamilton-Jacobi Bellman equations, J. Sci. Comput. 55 (3) (2013) 575–605.
  * [41] J. Garcke, A. Kröner, Suboptimal feedback control of pdes by solving hjb equations on adaptive sparse grids, J. Sci. Comput. 70 (2017) 1–28.
  * [42] W. Kang, L. C. Wilcox, Mitigating the curse of dimensionality: sparse grid characteristics method for optimal feedback control and HJB equations, Computational Optimization and Applications 68 (2017) 289–315.
  * [43] O. Bokanowski, X. Warin, A. Prost, Neural networks for first order hjb equations and application to front propagation with obstacle terms (2022).  arXiv:2210.04300.
  * [44] D. Kalise, K. Kunisch, Z. Rao, Hamilton-Jacobi-Bellman Equations: Numerical Methods and Applications in Optimal Control, Vol. 21 De Gruyter - Radon Series on Computational and Applied Mathematics, De Gruyter, Berlin, 2018\.
  * [45] K. Kunisch, D. Vásquez-Varas, Optimal polynomial feedback laws for finite horizon control problems, Comput. Math. with Appl. 148 (2023) 113–125.
  * [46] K. Kunisch, D. Vásquez-Varas, Smooth approximation of feedback laws for infinite horizon control problems with non-smooth value functions (2023).  arXiv:2312.11981.
  * [47] D. Onken, L. Nurbekyan, X. Li, S. W. Fung, S. Osher, L. Ruthotto, A neural network approach for high-dimensional optimal control applied to multiagent path finding, IEEE Trans. Control Syst. Technol. 31 (1) (2023) 235–251.
  * [48] K. Kunisch, D. Vásquez-Varas, D. Walter, Learning optimal feedback operators and their sparse polynomial approximations, J. Mach. Learn. Res. 24 (2023) 1–38.
  * [49] W. Sickel, T. Ullrich, Tensor products of Sobolev–Besov spaces and applications to approximation from the hyperbolic cross, J. Approx. Theory 161 (2) (2009) 748–786.
  * [50] W. H. Fleming, H. M. Soner, Controlled Markov Processes and Viscosity Solutions, second ed. Edition, Stochastic Modelling and Applied Probability, Springer, New York, 2006.
  * [51] P. Cannarsa, C. Sinestrari, Semiconcave Functions, Hamilton-Jacobi Equations, and Optimal Control, Birkhäuser Boston, 2004.

##  Appendix A Properties of $\mathcal{V}_{T}$

In this section we provide some important result concerning properties of
$\mathcal{V}_{T}$ needed throughout this work. The proofs can be obtained with
standard techniques and are therefore not presented here.

######  Lemma A1.

Let $T>0$ denote a finite time horizon, and let $y_{0}\in\omega$, $g_{1}\in
C^{2}(\Omega)$, $g_{2}\in C^{2}(\Omega)$, and $\delta>0$ be such that
$(T,u_{g_{1}},\delta,y_{0})$ and $(T,u_{g_{2}},\delta,y_{0})$ 3.1 holds. Then
we have for all $t\in[0,T]$

$|y(t;y_{0},u_{g_{1}})-y(t;y_{0},u_{g_{2}})|\leqslant\frac{|B|\left\lVert
u_{g_{1}}-u_{g_{2}}\right\rVert_{L^{\infty}(\Omega_{\delta};\mathbb{R}^{m})}}{a}\left(e^{ta}-1\right)$
(A.1)

and

$\begin{array}[]{ll}\displaystyle|\mathcal{V}_{T}(y_{0},u_{g_{1}})-\mathcal{V}_{T}(y_{0},u_{g_{2}})|\\\
\displaystyle\leqslant|B|\left(\frac{\left\lVert\nabla\ell\right\rVert_{L^{\infty}(\Omega_{\delta});\mathbb{R}^{d})}}{a}\left(\frac{e^{aT}-1}{a}-T\right)+\frac{T\beta}{2}\left\lVert
u_{g_{1}}+u_{g_{2}}\right\rVert_{L^{\infty}(\Omega_{\delta};\mathbb{R}^{d})}\right)\left\lVert
u_{g_{1}}-u_{g_{2}}\right\rVert_{L^{\infty}(\Omega_{\delta};\mathbb{R}^{m})},\end{array}$
(A.2)

where

$a=\left\lVert Df\right\rVert_{L^{\infty}(\Omega_{\delta};\mathbb{R}^{d\times
d})}+\frac{|B|^{2}}{\beta}\left\lVert\nabla^{2}g_{2}\right\rVert_{L^{\infty}(\Omega_{\delta};\mathbb{R}^{d\times
d})}.$

######  Lemma A2.

Let $T>0$ be a time horizon, an let $y_{0}\in\omega$, $u\in C^{1}(\Omega)$,
and $\delta>0$ be such that $(T,u,\delta,y_{0})$ 3.1 holds for all
$y_{0}\in\omega$. For arbitrary $y_{0}\in\omega$ let us define the adjoint
state $p(\cdot;y_{0},u)$ as the solution of the following adjoint equation

$-p^{\prime}(t)+\nabla\ell(y(t))+\beta Du(y(t))^{\top}\cdot
u(y(t))-(Df(y(t))+BDu(y(t)))^{\top}\cdot p(t)=0\mbox{ for all }t\in(0,T);\
p(T)=0,$ (A.3)

with $y(t)=y(t;y_{0},u)$. Then the function $\mathcal{V}_{T}$ its
differentiable with respect to each of its variables and we have

$D_{u}\mathcal{V}_{T}(y_{0},u)\cdot\nu=\int_{0}^{T}(\beta
u(y(t;y_{0},u))-B^{\top}p(t;y_{0},u))^{\top}\cdot\nu(y(t;y_{0},u))dt$ (A.4)

for all $\nu\in C^{1}(\Omega;\mathbb{R}^{m})$, and

$\nabla_{y}\mathcal{V}_{T}(y_{0},u)=-p(0;y_{0},u).$ (A.5)

Further, for all $t\in(0,T]$ we have

$\nabla_{y}\mathcal{V}_{T}(y(t),u)=-p(0;y(t),u)=-p(t;y_{0},u).$ (A.6)

######  Remark A.1.

Under the hypotheses of the previous lemma and assuming $u=u_{\phi}$ the
adjoint equation can be expressed in terms of $\phi$ as

$-p^{\prime}(t)+\nabla\ell(y(t))-Df^{\top}(y(t))p(t)+\frac{1}{\beta}\nabla^{2}\phi(y(t))BB^{\top}(\nabla\phi(y(t))+p(t))=0\mbox{
for all }t\in(0,T);\ p(T)=0.$ (A.7)

Additionally, using (A.6) we can express the derivative of
$\mathcal{V}_{T}(y_{0},u_{\phi})$ with respect to $\phi$ as

$D_{\phi}\left(\mathcal{V}_{T}(y_{0},u_{\phi})\right)\cdot\psi=\frac{1}{\beta}\int_{0}^{T}(\nabla\phi(y(t;y_{0},u_{\phi}))-\nabla_{y}\mathcal{V}_{T}(y(t;y_{0},u_{\phi})))^{\top}BB^{\top}\nabla\psi(y(t;y_{0},u_{\phi}))\,dt$
(A.8)

for $\psi\in C^{2}(\Omega)$.