%~Mouliné par MaN_auto v.0.27.3 2023-07-10 11:24:19
\documentclass[CRMECA, Unicode, XML, screen, Thematic, HideArticleType]{cedram}



\usepackage{graphicx}
\usepackage{mathtools}
\usepackage{booktabs}

\DeclareMathOperator{\supp}{supp}
\DeclareMathOperator*{\argmin}{arg\,min}
\DeclarePairedDelimiter{\abs}{\lvert}{\rvert}
\DeclarePairedDelimiter{\norm}{\lVert}{\rVert}






\newcommand{\eps}{\varepsilon}
\newcommand{\Om}{\Omega}

\newcommand{\N}{\mathbb{N}}
\newcommand{\R}{\mathbb{R}}
\newcommand{\de}{\mathrm{d}}
\newcommand{\arch}{\mathrm{arch}}
\newcommand{\net}{\mathrm{net}}
\newcommand{\spanl}{\mathrm{span}}
\newcommand{\col}{\mathrm{col}}
\newcommand{\Err}{\mathrm{Err}}
\DeclareMathOperator{\Errr}{Err}
\DeclareMathOperator{\intt}{int}
\DeclareMathOperator{\Symm}{Sym}
%\DeclareMathOperator{\arch}{arch}
\newcommand{\bfy}{\mathbf{y}}
\newcommand{\bfv}{\mathbf{v}}
\newcommand{\bfA}{\mathbf{A}}
\newcommand{\bff}{\mathbf{f}}
\newcommand{\bfg}{\mathbf{g}}
\newcommand{\bfp}{\mathbf{p}}
\newcommand{\bfQ}{\mathbf{Q}}
\newcommand{\bfU}{\mathbf{U}}
\newcommand{\bfY}{\mathbf{Y}}
\newcommand{\bfu}{\mathbf{u}}
\newcommand{\F}{\mathcal{F}}
\newcommand{\clB}{\mathcal{B}}
\newcommand{\clC}{\mathcal{C}}
\newcommand{\clY}{\mathcal{Y}}
\newcommand{\clO}{\mathcal{O}}
\newcommand{\clL}{\mathcal{L}}
\newcommand{\clV}{\mathcal{V}}
\newcommand{\clR}{\mathcal{R}}
\newcommand{\clJ}{\mathcal{J}}
\newcommand{\clI}{\mathcal{I}}
\newcommand{\clN}{\mathcal{N}}
\newcommand{\clZ}{\mathcal{Z}}
\newcommand{\M}{\mathcal{M}}
\newcommand{\clP}{\mathcal{P}}
\newcommand{\frkm}{\mathfrak{m}}

\newcommand{\Cbochnorm}[1]{\|#1\|_{\mathcal{C}}}
\newcommand{\wnorm}[1]{\|#1\|_{W_T}}

\newtheorem{assumption}{Assumption}
\renewcommand*{\theassumption}{\Alph{assumption}}

\newcommand{\tcb}{\textcolor{black}}
\newcommand{\black}{ \color{black} }
\newcommand{\bcdot}{ {\boldsymbol{\cdot}}}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\graphicspath{{./figures/}}

\newcommand*{\mk}{\mkern -1mu}
\newcommand*{\Mk}{\mkern -2mu}
\newcommand*{\mK}{\mkern 1mu}
\newcommand*{\MK}{\mkern 2mu}

%\hypersetup{urlcolor=purple, linkcolor=blue, citecolor=red}

\newcommand*{\relabel}{\renewcommand{\labelenumi}{(\theenumi)}}
\newcommand*{\romanenumi}{\renewcommand*{\theenumi}{\roman{enumi}}\relabel}
\newcommand*{\Romanenumi}{\renewcommand*{\theenumi}{\Roman{enumi}}\relabel}
\newcommand*{\alphenumi}{\renewcommand*{\theenumi}{\alph{enumi}}\relabel}
\newcommand*{\Alphenumi}{\renewcommand*{\theenumi}{\Alph{enumi}}\relabel}

\newcommand*{\Aenumi}{\renewcommand*{\theenumi}{\textbf{A.\arabic{enumi}}}\relabel}

\let\oldtilde\tilde
\renewcommand*{\tilde}[1]{\mathchoice{\widetilde{#1}}{\widetilde{#1}}{\oldtilde{#1}}{\oldtilde{#1}}}
\let\oldhat\hat
\renewcommand*{\hat}[1]{\mathchoice{\widehat{#1}}{\widehat{#1}}{\oldhat{#1}}{\oldhat{#1}}}
\let\oldforall\forall
\renewcommand*{\forall}{\mathrel{\oldforall}}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%\RubricEF{The scientific legacy of Roland Gloswinski}{L'héritage scientifique de Roland Gloswinski}
%
%
\title{Optimal feedback control of dynamical systems via value-function approximation}

\author{\firstname{Karl} \lastname{Kunisch}\IsCorresp}
\address{University of Graz, Institute of Mathematics and Scientific Computing, Heinrichstr. 36, A-8010 Graz, Austria}
\address{Johann Radon Institute for Computational and Applied Mathematics (RICAM), Austrian Academy of Sciences, Altenberger Stra\ss{}e 69, 4040 Linz, Austria}
\email[K. Kunisch]{karl.kunisch@uni-graz.at}

\author{\firstname{Daniel} \lastname{Walter}}
\address{Institut f\"ur Mathematik, Humboldt-Universit\"at~zu~Berlin, Rudower Chaussee 25, 10117 Berlin, Germany}
\email[D. Walter]{daniel.walter@hu-berlin.de}

\begin{abstract}
A self-learning approach for optimal feedback gains for finite-horizon nonlinear continuous time control systems is proposed and analysed. It relies on parameter dependent approximations to the optimal value function obtained from a family of universal approximators. The cost functional for the training of an approximate optimal feedback law incorporates two main features. First, it contains the average over the objective functional values of the parametrized feedback control for an ensemble of initial values. Second, it is adapted to exploit the relationship between the maximum principle and dynamic programming. Based on universal approximation properties, existence, convergence and first order optimality conditions for optimal neural network feedback controllers are proved.
\end{abstract}

\keywords{optimal feedback control, neural networks, Hamilton--Jacobi--Bellman equation, self-learning, reinforcement learning}

\subjclass{49J15, 49N35, 68Q32, 93B52, 93D15}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\dateposted{2023-07-20}
\begin{document}
\maketitle


\section{Introduction} \label{sec:Introduction}

In this paper we focus on optimal feedback control for problems of the form
\begin{equation}\label{def:refproblem}
\tag{$P$}
\begin{cases}
\quad &\inf_{y,u} J(y,u) \coloneqq \frac{1}{2}\int_{0}^T \left(|Q_1(y(t)-y_d(t))|^2+ \beta |u(t)|^2 \right)~\de t +\frac{1}{2}\left|Q_2\left(y(T)-y^T_d\right)\right|^2\\
&s.t. \quad \dot{y}= {f}(y)+g(y)u, \quad y(0)=y_0, \text{ and } u \in L^2\left(0,T; \R^m\right),
\end{cases}
\end{equation}
with nonlinear dynamics described by $f:[0,T]\times \R^n \to \R^n$. The system can be influenced by choosing a control input~$u$ which enters through a control operator $g:\R^n\to \R^{n \times m}$. We assess the perfomance of a given control by its objective functional value which comprises the (weighted) distance between the associated state trajectory~$y$ and a given desired state~$y_d$ as well as the norm of the control for some cost parameter~$\beta>0$. The weighting matrices $Q_i$, for $i=1,2$, are assumed to be symmetric positive semi-definite. Searching for an optimal control $u^*$ in feedback form requires to find a function $F^*:[0,T]\times \R^n \to \R^m$ such that
\[
u^*(t)= F^*\left(t,y^*(t)\right), \text{ for } t\in (0,T).
\]
Here $(u^*,y^*)$ denotes an optimal control-trajectory pair associated to~\eqref{def:refproblem}. Under appropriate conditions, see e.g.~\cite{FS06}, the feedback mapping can be expressed as
\begin{equation}\label{eq:intro1}
F^*(t,y) = -\frac{1}{\beta}g^\top(t,y) \partial_y V^*(t,y),
\end{equation}
where $V^*$ stands for the value function associate to~\eqref{def:refproblem}, i.e. for $(T_0,y_0)\in [0,T]\times \R^n$:
\begin{align*}
V^*(T_0,y_0) &= \min_{y,u} J_{T_0}(y,u), \text{ subject to } \dot{y}= {f}(y)+g(y)u, \quad y(T_0)=y_0,
\\
\intertext{and}
J_{T_0}(y,u) &=\frac{1}{2}\int_{T_0}^T \left (\left|Q_1\left(y(t)-y_d(t)\right)\right|^2+ \beta |u(t)|^2 \right)~\de t +\frac{1}{2}\left|Q_2\left(y(T)-y^T_d\right)\right|^2.
\end{align*}
The value function $V^*$ satisfies a Hamilton--Jacobi--Bellman (HJB) equation which is a time-dependent first order hyperbolic equation of spatial dimension $n$. Numerical realisations, therefore, are plagued by the curse of dimensionality. Indeed a direct solution of the HJB equation already becomes computationally prohibitive for moderate dimensions~$n$.


Therefore, for practical realization, the interest in alternative techniques arises. In many situations of practical relevance researches have relied on linear approximations to the nonlinear dynamical system and have treated the resulting linear-quadratic problem by Riccati techniques. Much research has concentrated on validating this approach locally around a reference trajectory. Globally such a strategy may fail, see for instance~\cite{dkk19, KW2020}.

In this paper we follow an approach, possibly first proposed in~\cite{KW2020}, circumventing the construction of the value function on the basis of solving the HJB equation. Rather the feedback mapping is constructed by an unsupervised self-learning technique. In practice, this requires the approximation of~$V^*$ by a family of functions~$V_\theta$ which are parametrized by a finite dimensional vector~$\theta$ and satisfy a uniform approximation property. Possible families of universal approximators include, e.g., neural networks or piecewise polynomial approximations. Subsequently, in view of~\eqref{eq:intro1}, we introduce the corresponding feedback law
\begin{equation}\label{eq:intro2}
F_\theta(t,y) = -\frac{1}{\beta}g^\top(y) \partial_y V_\theta(t,y), \text{ for } (t,y)\in [0,\infty) \times \R^n,
\end{equation}
as approximation to $F^*$. An ``optimal'' parametrized feedback law is then determined by a variant of the following self-learning, structure preserving, variational problem:
\begin{multline}\label{eq:intro3}
\min_{\theta}\; J(y,\F_{\theta}(y))
+\frac{1}{2}\int^T_0 {\gamma_1} \left|V_{\theta}(t,y(t))-J_t\left(y,F_{\theta}(\cdot,y)\right)\right|^2 + {\gamma_2}\left|\partial_y V_{\theta}(t,y(t))-p(t)\right|^2~\de t +\frac{\gamma_\eps}{2} {|\theta|^2}\\
\begin{aligned}
\text{s.t. } \quad \dot{y} &= {f}(y)+g(y) F_{\theta}(y), \quad y(0)=y_0, \quad p(T)= Q^\top_2 Q_2\left(y(T)-y^T_d\right)\\[1.5ex]
- \dot p &= {f}(y)^\top p+\left| D {g}(y)^\top F_\theta (y)\right| p+Q_1^\top Q_1(y-y_d).
\end{aligned}
\end{multline}
In this problem, minimization with respect to $u$ is replaced by minimizing with respect to the parameters $\theta$ which characterize $V_\theta$ and~$F_\theta$. The cost functional of problem~\eqref{eq:intro3} consists of four parts: The first term represents the objective functional of~\eqref{def:refproblem} where the control~$u$ is replaced by the closed loop expression~$F_\theta(y)$. The next two terms realize the fact that $V_\theta$ is constructed as approximation to the value function associated to~\eqref{def:refproblem} and exploit the well-known property that, under certain conditions, the gradient of the value function coincides with the solution of a suitable adjoint equation, see e.g.~\cite[p.~21]{FS06}. The final term penalizes the norm of the structural parameters. We point out that $V_\theta$ and $F_\theta$ are learned along the orbit $\clO= \{y(t;y_0): t\in (0,\infty)\}$ within the state space $\R^n$. To accommodate the case that one trajectory does not provide enough information, we propose to involve an ensemble of orbits departing from a set $Y_0$ of initial conditions, and to reformulate problem~\eqref{eq:intro3} accordingly. This will be done in Section~\ref{sec:learnfeedback} below. While we focus on linear-quadratic objective functionals, the derived results can be readily generalized if suitable coercivity and differentiability properties are assumed for this functional.

In our earlier work on learning a feedback function~\cite{KW2020}, we considered infinite horizon optimal control problems. In that case, the time-dependent HJB equation results in a stationary one. There we had not yet incorporated the structure preserving terms involving $V_\theta$ and $\partial_y V_\theta$ into the cost. Moreover we directly constructed an approximation $F_\theta$ to the vector valued function $F^*$, rather than approximating the scalar valued function $V^*$ and subsequently using~\eqref{eq:intro2}. In the present paper we provide the theoretical foundations for the learning based technique that we propose to construct an approximation to the optimal feedback function for~\eqref{def:refproblem}. Recently in~\cite{onlfor2020} a variant of the approach as in~\cite{KW2020} was used for interesting numerical investigations to construct optimal feedback functions for finite horizon multi-agent optimal control problems.

Let us very briefly mention some of the vast literature on solving the HJB equations. Semi-Lagrangian schemes and finite difference methods have been deeply investigated to directly solve HJB equations directly, see e.g.~\cite{bggk13, ff16, kkr18}. Significant progress was made in solving high dimensional HJB equations by the of use policy iterations combined with tensor calculus techniques,\linebreak\cite{dkk19, kk18, foss2020}. The use of Hopf formulas was proposed in e.g.~\cite{lr86, cloy18}. Interpolation techniques, utilizing ensembles of open loop solutions have been analyzed in the works of~\cite{akk2020, ngk19}, for example. Finally we mention that optimal feedback control is intimately related to reinforcement learning, see e.g. the monograph~\cite{bert19}, and also the survey articles~\cite{LV09, recht18, vls14}.



The manuscript is structured as follows. Some pertinent notation is gathered in Section~\ref{sect2}. In Section~\ref{sec:o} concepts of optimal feedback control, semi-global with respect to the initial condition $y_0$, are gathered. Section~\ref{sec:learnfeedback} is devoted to describing the learning technique that we propose to approximate the optimal feedback function. In Section~\ref{sec:examples} the required assumptions on approximating subspaces are checked for a class of neural networks and a class of piecewise polynomials. Existence of solutions to the approximating learning problems is proved in Section~\ref{sec:existence}. Their convergence is analyzed in Section~\ref{sec:convergence}. The case of learning from finitely many orbits is the focus of Section~\ref{sect8}. Section~\ref{sec:numericalexampl} provides an example illustrating the numerical feasibility of the proposed method. We do not aim for sophistication in this respect. The Appendix~\ref{app} details the proofs of several necessary technical results.



\section{Notation}\label{sect2}

For~$I:=(0,T)$, with $T>0$, we define $ W_T= \{\,y \in L^2(I; \R^n)\;|\;\dot{y}\in L^2(I; \R^n)\,\}, $ where the temporal derivative is understood in the distributional sense. We equip~$W_T$ with the norm induced by the inner product
\begin{align*}
(y_1,y_2)_{W_T}=(\dot{y}_1,\dot{y}_2)_{L^2(I;\R^n)}+(y_1,y_2)_{L^2(I;\R^n)} \quad \text{for } y_1, y_2 \in W_T,
\end{align*}
making it a Hilbert space. We recall that $W_T$ embeds continuously into $C(\bar I; \R^n)$. For a compact metric space~$X$ we denote the space of continuous functions between~$X$ and~$Y$ by~$\clC (X;Y)$ which we endow with $
\|\varphi\|_{\clC (X;Y)}= \max_{x \in X} \|\varphi(x)\|_Y$ as norm. By $Y_0$ we denote a compact set of initial conditions in $\R^n$. When arising as index, the space $\clC (Y_0;W_T)$ will frequently be abbreviated by $\clC $. The space~$\clC ^1(X;Y)$ of continuously differentiable functions is defined analogously. Open balls of radius $\eps$ in a Banach space $X$ with center $x$ will be denoted by $B_{\eps}(x)$. The space of bounded linear operators between Banach spaces~$X$ and~$Y$, endowed with the canonical norm, is denoted by~$\clB (X,Y)$. We further abbreviate~$\clB (X):=\clB (X,X)$.

\section{Semi-global optimal feedback control}\label{sec:o}
Consider the controlled nonlinear dynamical system of the form
\begin{align}\label{eq:openloop}
\dot{y}= \bff(y)+\bfg(y) u \quad \text{in}~L^2\left(I;\R^n\right), \quad y(0)=y_0,
\end{align}
described by Nemitsky operators
\begin{equation}\label{eq:spaces}
\begin{split}
\bff\colon W_T &\to L^2\left(I;\R^n\right), \mkern140mu \bff(y)(t)= f(t,y(t)) \\[1.3ex]
\bfg \colon W_T &\to \clL\left(L^2\left(I;\R^m\right);L^2\left(I;\R^n\right)\right), \quad \bfg(y)(t)= g(t,y(t))
\end{split}
\end{equation}
for~a.e.~$t\in I$,~$f \colon I \times \R^n \to \R^m$ and~$g \colon I \times \R^n \to \R^{n \times m}$. The smoothness requirements on~$f$ and~$g$ will be detailed in Assumption~\ref{ass:feedbacklaw} below. Our aim is to choose a control input~$u^*\in L^2(I;\R^m)$ which keeps the associated solution~$y^* \in W_T$ close to a known reference trajectory~$y_d$, while keeping the control effort small. This is formulated as the constrained minimization problem
\begin{equation}\label{def:openloopproblem}
\tag{$P_{y_0}$}
\begin{cases}
\quad \inf\limits_{y\,\in\,W_T,\, u\,\in\,L^2(I; \R^m)} J(y,u) \\
s.t. \quad \dot{y}= \bff(y)+\bfg(y) u, \quad y(0)=y_0,
\end{cases}
\end{equation}
where
\[
J(y,u)=\frac{1}{2}\int_{I} \left (\left|Q_1\left(y(t)-y_d(t)\right)\right|^2+ \beta |u(t)|^2 \right)~\de t+\frac{1}{2}\left|Q_2\left(y(T)-y^T_d\right)\right|^2,
\]
which incorporates the weighted misfit between the trajectory~$y$ within the time horizon $I=(0,T)$ and at the terminal time to desired states~$y_d\in L^2(I;\R^n)$ and $y_d^T \in \R^n$, as well as the norm of the control~$u$. While this~\emph{open loop} optimal control problem captures well the objective formulated above, it comes with several disadvantages. First, its solution is a function of time only, and does not include the current state $y(t)$. This makes the open loop approach susceptible to possible perturbations in the dynamical system. Second, determining the control action for a new initial condition requires to solve~\eqref{def:openloopproblem} from the start.

The aforementioned limitations of open loop optimal controls motivate the study of~\emph{semi-global optimal feedback control} approaches to~\eqref{def:openloopproblem}. More precisely, given a compact set~$Y_0 \subset \R^n $, we look for a feedback function~$F^* \colon I \times \R^n \to\R^m$ which induces a Nemitsky operator
\begin{align*}
\F^* \colon W_T \to {L^2}\left(I; \R^m\right), \quad \F^*(y)(t)=F^*(t,y(t)) \quad \text{for a.e.}~t \in I,
\end{align*}
such that for every~$y_0 \in Y_0$ the~\emph{closed loop system}
\begin{equation}\label{eq:cloloop}
\dot{y}= \bff(y)+ \bfg(y) \F^*(y), \quad y(0)=y_0,
\end{equation}
admits a unique solution~$y^*(y_0) \in W_T$ and~$(y^*(y_0), \F^*(y^*(y_0)))$ is a minimizing pair of~\eqref{def:openloopproblem}.


The determination of an optimal feedback function usually rests on the computation of the value function to\eqref{def:openloopproblem} which is defined as
\begin{align}\label{def:valuefunc}
V^*(T_0,y_0):= \min_{\substack{y\,\in\,H^1\left(T_0,T;\R^n\right), \\ u \in L^2\left(T_0,T; \R^m\right)}} J_{T_0}(y,u) \quad s.t. \quad \dot{y}= \bff(y)+ \bfg(y) u, \quad y(t_0)=y_0,
\end{align}
where~$(T_0,y_0)\in I \times \R^n$, and $J_{T_0}(y,u)$ is defined as
\begin{align*}
J_{T_0} (y,u)= \frac{1}{2}\int^T_{T_0} \left(\left|Q_1(y(t)-y_d(t))\right|^2+ \beta |u(t)|^2 \right)~\de t+\frac{1}{2}\left|Q_2(y(T)-y_d(T))\right|^2.
\end{align*}
By construction $V^*$ satisfies the final time boundary condition
\begin{align*}
V^*(T,y_0)=\frac{1}{2} \left|Q_2(y_0-y_d(T))\right|^2 \quad \forall y_0 \in \R^n.
\end{align*}
If~$V^*$ is continuously differentiable in a neighborhood of some~$(t,y_0) \in I \times \R^n$ then it solves the instationary~\emph{Hamilton--Jacobi--Bellman (HJB) equation}
\begin{align}\label{eq:HJB}
\partial_t V^*(t,y_0) + \left(f(y_0), \partial_y V^*(t,y_0)\right)_{\R^n}- \frac{1}{2\beta}\left|g(t,y_0)^\top\partial_y V^*(t,y_0)\right|^2+ \frac{1}{2} \left|Q_1 (y_0-y_d(t))\right|^2=0
\end{align}
in the classical sense there, see e.g.~\cite{ff14,FS06}. Here~$\partial_t V^*$ denotes the partial derivative of the value function with respect to~$t$ and~$\partial_y V^*$ is the gradient of~$V^*$ with respect to the~$y$-variable. An optimal control for~\eqref{def:openloopproblem} in feedback form is then given by $u^*=-\frac{1}{\beta} \bfg(y^*)^\top\partial_y \clV^*(y^*)$ where~$\partial_y \clV^*(y^*)(t)=\partial_y V^*(t,y^*(t))$ for every~$t\in I$, and~$y^*=y^*(y_0) \in W_T$ solves the closed loop system
\begin{align*}
\dot{y}=\bff(y)- \frac{1}{\beta}\bfg(y) \bfg(y)^\top \partial_y \clV^*(y), \quad y(0)=y_0.
\end{align*}
Thus
\begin{align*}
\left(y^*(y_0), -\frac{1}{\beta}\bfg\left(y^*(y_0)\right)^\top \partial_y \clV^*\left(y^*(y_0)\right)\right) \in \argmin~\eqref{def:openloopproblem}
\end{align*}
and the function
\begin{align*}
F^*(\cdot,\cdot)=-\frac{1}{\beta} g(\cdot,\cdot)^\top \partial_y V^*(\cdot,\cdot)
\end{align*}
is an optimal feedback law.

Realizing the optimal feedback in this way requires a solution to~\eqref{eq:HJB} which is a partial differential equation on~$\R^n$. This can be extremely challenging or even impossible depending on the dimension $n$ and the computational facilities at hand. Similarly to our previous manuscript~\cite{KW2020}, we take a different approach by formulating minimization problem over a suitable set of feedback functions involving the closed loop system as a constraint. This relates to a learning problem, within which the feedback functions are trained to achieve optimal stabilization. This makes the problem computationally amenable.


The procedure just described will be formalized in the following section. Here we first summarize the assumptions on the nonlinear dynamical system that we refer to throughout the paper.


\begin{enonce}{Assumption}\label{ass:feedbacklaw}
\leavevmode
\begin{enumerate}\Aenumi
%[(\textbf{A.1})]
\item\label{A.1}The functions~$f \colon I \times \R^n \to \R^n$ and~$g \colon I \times \R^{n} \to \R^{n \times m}$ are twice continuously differentiable. Their Jacobians and Hessians with respect to the second variable, denoted by~$D_{y}f, D_{yy}f$, and~$D_{y}g, D_{yy}g$, respectively, are Lipschitz continuous on compact sets, uniformly for $t\in I$.
%[(\textbf{A.2})]
\item\label{A.2} There exists a constant~$M_{Y_0}>0$ such that the value function~$V(\cdot,\cdot)\colon I \times \R^n \to \R$ for~\eqref{def:openloopproblem} is twice continuously differentiable on~$I \times \bar{B}_{2\widehat{M}}(0)$ with Lipschitz continuous gradient and Hessian (w.r.t. $y$ uniformly in~$t \in I$) where
\begin{align}\label{eq:aux1}
\widehat M = M_{Y_0} \, \|\imath\|_{\clB\left(W_T,\, \clC (I;\R^n)\right)},
\end{align}
and $\imath$ denotes the embedding of $W_T$ into $\clC (I;\R^n))$.
\end{enumerate}
\end{enonce}


As a consequence of~\eqref{A.1}, the Nemitsky operators~$\bff,~\bfg$ are at least two times continuously differentiable with domains and ranges as defined in~\eqref{eq:spaces}. Their derivatives, denoted by $D \bff(y) \in{\clB}(W_T,L^2(I;\R^n)),$ and $D \bfg(y) \in \clB (W_T; \clB ({L^2(I;\R^m); L^2(I;\R^n)})),$ are the Nemitsky operators induced by~$D_y f$ and~$D_y g$. Moreover~$ \bff, D\bff, \bfg,D \bfg$ are Lipschitz continuous and bounded, on bounded subsets of $L^\infty(I;\R^n)$, and thus in particular on~$\clY _{ad}\subset W_T$, where
\begin{align}\label{eq:kk13}
\clY _{ad}:= \left \{\,y \in W_T\;|\;\wnorm{y}\leq 2 M_{Y_0}\,\right\}.
\end{align}
Finally $D\bff^\top(y) \in{\clB}(W_T,L^2(I;\R^n))$ denotes the Nemitsky operator associated to $D_y f^\top$.

Analogously, due to~\eqref{A.2}, $V^*$ induces a twice Lipschitz continuously Fréchet differentiable Nemitsky operator~$\clV^*: \clY _{ad}\subset W_T \to L^{\tcb{2}}(I)$. Moreover~$\clV^*$ and its first derivative~$D\clV^* $ are weak-to-strong continuous. Define the Nemitsky operator
\begin{align}\label{def:optfeddnemitsk}
\F^* \colon \clY _{ad}\to L^2\left(I;\R^n\right), \tcb \quad \F^*(y)=-\frac{1}{\beta} \bfg(y)^\top \partial_y \clV^*(y),
\end{align}
where~$\partial_y \clV^*$ is the Nemitsky operator induced by the gradient~$\partial_y V^*=D_y V (\cdot,\cdot)^\top$. Note also that $\F^* \in C^1(W_T; (L^2(I;\R^m);L^2(I;\R^n)))$. We further assume the following:
\begin{enumerate}\Aenumi
\setcounter{enumi}{2}
\item\label{A.3} For every~$y_0 \in Y_0$ there exists a unique function~$y=\bfy^*(y_0)\in W_T$ satisfying
\begin{align*}
\dot{y}= \bff(y)+ \bfg(y) \F^*(y), \quad y(0)=y_0, \quad \wnorm{y}\leq M_{Y_0}.
\end{align*}
Moreover we have
\begin{align*}
\left(y^*(y_0), \F^*\left(y^*(y_0)\right)\right) \in \argmin~\eqref{def:openloopproblem} \quad \forall y_0 \in Y_0.
\end{align*}
\end{enumerate}
When referring to Assumption A we mean~\eqref{A.1}-\eqref{A.3}. We emphasize that the constant $M$ appearing in~\eqref{A.2} and~\eqref{A.3} is assumed to be same. Note further that as a consequence of~\eqref{A.3} problem~\eqref{def:openloopproblem} admits a solution for each $y_0 \in Y_0$, with the optimal control given by $u^*= \F^*(y^*(y_0)))$.

\begin{rema}\label{rem3}
Using~\eqref{A.1}, \eqref{A.3} as well as the implicit function theorem it can be readily be verified that the mapping~$\bfy^*\colon Y_0 \to W_T$ from~\eqref{A.3} is continuously differentiable. Given~$\delta y_0 \in \R^n$ the directional derivative~$\delta y\coloneqq \partial\bfy^*(y_0)(\delta y_0)$ of~$\bfy^*$ at~$y_0 \in Y_0$ in direction $\delta y_0$ satisfies the linearized ODE system
\begin{align*}
\dot{\delta y}= D \bff\left(\bfy^*(y_0)\right)\delta y+ \left\lbrack D\bfg\left(\bfy^*(y_0)\right)\delta y\right\rbrack \F^*\left(\bfy^*(y_0)\right) +\bfg\left(\bfy^*(y_0)\right)D\F^*\left(\bfy^*(y_0)\right)\delta y,\, \delta y(0)= \delta y_0.
\end{align*}
Here $D \bfg$ is induced by $D_yg$ which is given by
\begin{align*}
\left \lbrack D_y g(t,y) \delta y \right\rbrack_{ij}= \left(\sum^n_{k=1} \partial_k g_{ij}(t,y) \delta y_k \right) \quad \forall \delta y \in \R^n,
\end{align*}
where~$g(y)=(g_{ij})$ and~``$\partial_k$'' denotes the partial derivative w.r.t to the~$k^{\rm th}$ component of~$y$. The transposed~$D\bfg(y)^\top$, which will arise in the adjoint equation below, is induced by the tensor~$D_y g(t,\cdot)^\top=(D_y g(t,\cdot)_{kji})\in \R^{n\times n\times m}$, with $t\in I$. In particular, we readily verify that~$D \bfg(\cdot)^\top \in \clB (L^2(I;\R^m); \clB (W_T; L^2(I;\R^n)))$.
\end{rema}


To end this section we collect structural information on the relation between the adjoined state, denoted by ${p}$ below, the optima value function~$V^*$, and the induced optimal feedback law~$\F^*$.

\begin{prop}\label{prop:structure}
Let Assumption~\ref{ass:feedbacklaw} hold. Then there exists a unique continuous mapping~$\bfp^* \colon Y_0 \to W_T$ such that for each~$y_0 \in Y_0$ the tuple~$(y,p)=(\bfy^*(y_0),\bfp^*(y_0))$ satisfies
\begin{align}
\frac{d}{dt} y&=\bff(y)+\bfg(y)\F^*(y),~y(0)=y_0, \label{eq:stateprop}
\\
-\frac{d}{dt} p&= D\bff(y)^\top p+\left\lbrack D \bfg(y)^\top\F^*(y)\right\rbrack p+Q_1^\top Q_1(y-y_d),~p(T)= Q^\top_2 Q_2\left(y(T)-y^T_d\right), \label{eq:adjointprop}
\\
\F^*(y)&=- \frac{1}{\beta} \bfg(y)^\top p. \label{eq:gradienteqprop}
\end{align}
Moreover we have
\begin{align}\label{eq:dynamicalprog}
V^*(t,y(t))=J_t\left(y(t),F^*(t,y(t))\right),~p(t)= \partial_y V(t,y(t)) \quad \forall t\in [0,T].
\end{align}
\end{prop}

\begin{proof}[Proof of Proposition~\ref{prop:structure}]
By~\eqref{A.3} problem~\eqref{def:openloopproblem} admits a solution for each $y_0 \in Y_0$. Then~\eqref{A.1}-\eqref{A.2} guarantee that~\eqref{eq:adjointprop}, with $y=y(y_0)\in W_T$ the state component of a solution to~\eqref{def:openloopproblem}, admits a unique solution $p$ in~$W_T$ which continuously depends on~$y\in W_T$. Moreover~\eqref{eq:stateprop} -~\eqref{eq:gradienteqprop} represent the first order necessary optimality condition for~\eqref{def:openloopproblem} with the optimal control $u(t)= \F^*(y(t))$. Since $\bfy^*\colon Y_0 \to W_T$ is continuous as mentioned in Remark~\ref{rem3} and the solution to~\eqref{eq:adjointprop} depends continuously on $y\in W_T$, the claimed continuity $\bfp^* \colon Y_0 \to W_T$ follows. Equation~\eqref{eq:dynamicalprog} is a direct consequence of the dynamic programming principle, and~\eqref{A.3}.
\end{proof}



\section{Optimal feedback control by value function approximation}\label{sec:learnfeedback}
This section is devoted to introducing a family of computationally tractable minimization problems from which we will ``learn'' approximations of optimal feedback laws. Our approach rests on two main pillars. First, given~$\eps>0$, we consider a family of functions $V^\eps_{\theta} \in \clC (I\times \R^n) $ which are finitely parametrized by~$\theta \in \clR_\eps \simeq \R^{N_\eps}$,~$N_\eps \in \N$. These serve as ``discrete'' approximations of the optimal value function~$V^*$. The following a priori estimate is assumed, for some fixed $\eps_0 >0$:

\begin{enonce}{Assumption} \label{ass:approxsmoothness}
For every~$0 <\eps \leq \eps_0$ there holds~$V_{\cdot}^\eps \in\clC ^4(\clR_\eps \times \R \times \R^n)$ and~$V^\eps_{\theta}(T,y_0)=\frac{1}{2} |Q_2(y_0-y_d(T))|^2$ for every~$y_0 \in \R^n$ and~$\theta \in \clR_\eps$. Moreover there exists~$\theta_\eps \in \clR_{\eps}$ with
\begin{align}\label{eq:approxcapa}
\max_{\substack{t\,\in\,I, \\ |y|\,\leq\,2 \widehat{M} }} \left|V^\eps_{\theta_\eps}(t,y)-V^*(t,y)\right|+\left|\partial_y \left(V^\eps_{\theta_\eps}(t,y)-V^*(t,y)\right)\right| +\left\|\partial_{yy}\left(V^\eps_{\theta_\eps}(t,y)-V^*(t,y)\right)\right\| \leq c \eps
\end{align}
for some~$c>0$ independent of~$\eps\in (0,\eps_0]$.
\end{enonce}

Now recall from~\eqref{def:optfeddnemitsk} that the optimal feedback law~$\F^*$ is the superposition operator induced by~$F^*(t,y)=-(1/\beta)g(t,y)^\top \partial_y V^*(t,y) $. With the aim of preserving the dependence of the feedback law on the value function in our approximation, we define a set of parametrized feedback laws~$\F^\eps_\theta$ associated to~$V^\eps_\theta$,~$\theta \in \clR_\eps$, by
\begin{align*}
\F^\eps_\theta(y)(t)= F^\eps_\theta(t,y(t))= -\frac{1}{\beta} g(t,y(t)) \partial_y V^\eps_\theta(t,y(t))
\end{align*}
for all~$y\in W_T$,~$t \in \bar{I}$ and~$\theta \in \clR_\eps$. A first approach to obtain an optimal feedback law in the form~$\F^\eps_\theta$ can then be found by replacing the open loop control~$u$ in~\eqref{def:openloopproblem} by the closed loop expression~$\F^\eps_\theta(y)$ and minimizing for~$\theta \in \clR_\eps$:
\begin{align}\label{eq:aux7}
\min_{y\,\in\,W_T, \theta\,\in\,\clR_\eps} J\left(y, F^\eps_\theta(y)\right)+ \frac{\gamma_\eps}{2} \|\theta\|^2_{\clR_\eps} \quad s.t. \quad \dot{y}= \bff(y)+\bfg(y) \F^\eps_\theta(y), y(0)=y_0,
\end{align}
where~$\|\cdot\|_{\clR_\eps}$ denotes a Hilbert space norm on~$\clR_\eps$,~$\gamma_\eps > 0$ and~$y_0 \in Y_0$ is fixed. This represents the goal of finding a feedback law~$\F^\eps_\theta$ together with a trajectory~$y\in W_T$ which satisfy~$(y, \F^\eps_\theta(y))\in \argmin~\eqref{def:openloopproblem} $. However, this approach falls short in several aspects. First, we cannot hope to recover a solution of the semiglobal optimal feedback control problem for all $y_0\in Y_0$, since the minimization in~\eqref{eq:aux7} is associate to a single initial condition only. Secondly it misses to impose properties that would guide $\F^\eps_\theta(y)$ to be close to $\clV^*$, and it does not exploit the relation between the adjoint state~$p$, see~\eqref{eq:adjointprop}, and the gradient of the value function~$\partial_y \clV^*$. Incorporating this information into the problem can, potentially, lead to improved learning results and improved parameterized feedback laws which behave similarly to~$\F^*$. These considerations lead to the second pillar of our approach, namely a succinct choice of the cost for the learning problem. For this purpose we use all of $Y_0$ as ``learning set'' for initial conditions. It is endowed with the normalized Lebesgue measure~$\clL$. Moreover we define the augmented objective
\begin{equation}\label{eq:aux2}
J_\eps(y,p,\theta)= J\left(y,\F^\eps_{\theta}(y)\right)+\int^T_0 \frac{\gamma_1}{2} \left|V^\eps_{\theta}(t,y(t))-J_t\left(y,\F^\eps_{\theta}(y)\right)\right|^2 + \frac{\gamma_2}{2}\left|\partial_y V^\eps_{\theta}(t,y(t))-p(t)\right|^2~\de t
\end{equation}
for penalty parameters~$\gamma_1, \gamma_2 \geq 0$. The arguments in $J_t$ are the restriction of the solution $y$ to the equation in~\eqref{eq:aux7} and the feedback $\F^\eps_{\theta}(y)$ to $[t,T]$. The additional terms in this new objective functional penalize the violation of the cost and its gradient by means of the approximation based on $V_\theta^\eps$, i.e. they penalize the differences between $J_t(y,\F^\eps_\theta(y))$ and $V^\eps_\theta(t,y(t)$, as well as $p(t)$ and~$\partial_y V^\eps_\theta(t,y(t)$.


Given a strictly positive weight function~$\omega \in L^\infty(Y_0);~0<c \leq \omega$ a.e., we thus propose to find a feedback law~$\F^\eps_\theta$ by solving the ensemble control problem
\begin{align}\label{def:approxfeedprop}
\min_{\substack{\bfy\,\in\,\bfy_{ad},\\ \F^\eps_\theta(\bfy)\,\in\,\bfU_{ad},\\
\bfp\,\in\,\clC (Y_0;W_T) \\
\theta\,\in\,\clR_{\eps}}}\;\clJ_\eps(\bfy,\bfp,\theta)\coloneqq \int_{Y_0} \omega(y_0)\, J_\eps\left(\bfy(y_0),\bfp(y_0),\theta\right)\,\de \clL(y_0)+ \frac{\gamma_\eps}{2} \|\theta\|^2_{\clR_\eps} \tag{$\clP_\eps$}
\end{align}
subject to the system of closed loop state~\emph{and} adjoint equations
\begin{align}
\dot{\bfy}(y_0) &=\bff(\bfy(y_0))+\bfg(\bfy(y_0))\F^\eps_\theta(\bfy(y_0)) \label{eq:statepropapprox}
\\
-\dot{\bfp}(y_0) &= D\bff(\bfy(y_0))^\top\bfp(y_0)+ \left\lbrack D\bfg(\bfy(y_0))^\top\F^\eps_\theta(\bfy(y_0))\right\rbrack \bfp(y_0)+\bfQ_1^\top \bfQ_1\left(\bfy(y_0)-y_d\right) \label{eq:adjointpropapprox}
\\
\bfy(y_0)(0) &=y_0, \bfp(y_0)(T)= Q^\top_2 Q_2\left(\bfy(y_0)(T)-y_d^T\right), \bfy(y_0)\in \clY _{ad} \label{eq:constraintprop}
\end{align}
for~$\clL$-a.e.~$y_0\in Y_0$. Above~$\bfy_{ad} \subset \clC (Y_0;W_T)$ and~$\bfU_{ad} \subset L^2(Y_0;L^2(I;\R^m))$ denote the admissible sets of ensemble state trajectories and admissible controls. They will be specified in section~\ref{sec:existence}.



\section{Examples} \label{sec:examples}
In this section we discuss two particular examples for the parameterized mappings~$V^\eps$: deep residual networks and piecewise polynomial functions of sufficiently high degree.


\subsection{Residual networks} \label{subsec:residual}
To explain the approximation of the value function by residual neural networks, we first fix some notation. Let~$L_\eps\in \N$,~$L_\eps \geq 2$, as well as~$N^{\eps}_{i}\in
\N$, $i=1,\,\dots,\,L_\eps -1$ be given. We set~$N^{\eps}_0=n+1$ and~$N^{\eps}_L=1$. Furthermore define
\begin{align*}
\clR_\eps= \bigtimes^{L_\eps-1}_{i=1} \left (\R^{N^\eps_{i}
\times N^\eps_{i-1}} \times \R^{N^\eps_{i} \times N^\eps_{i-1}} \times
\R^{N^\eps_i} \right) \times \R^{N^\eps_{L} \times N^\eps_{L-1}}.
\end{align*}
The space~$\clR_\eps$ is uniquely determined by its~\emph{architecture}
\begin{align*}
\arch(\clR_\eps)=\left(N^\eps_0, N^\eps_1,\,\dots,\, N^\eps_L\right)\in \N^{L_\eps+1}.
\end{align*}
A set of parameters~$\theta \in \clR_\eps$ given by
\begin{align*}
\theta=\left(W_{11}, W_{12}, b_1,\,\dots,\, W_{{L_\eps}}
\right)
\end{align*}
is called a~\emph{neural network} with~$L_\eps$~\emph{layers}. Moreover let~$\sigma \in \clC ^4(\R)$ be given and assume that~$\sigma$ is not a polynomial. The function
\begin{equation}\label{eq:shift}
V^\eps_{\theta}(t,y)= \frac{1}{2} \left|Q_2\left(y-y_d(T)\right)\right|^2 + f^{\sigma}_{L_\eps, \theta} \circ f^{\sigma}_{L_\eps-1, \theta} \circ \cdots \circ f^{\sigma}_{1,
\theta}((t,y))-f^{\sigma}_{L_\eps, \theta} \circ f^{\sigma}_{L_\eps-1,
\theta} \circ \cdots \circ f^{\sigma}_{1, \theta}((T,y))
\end{equation}
for~$(t,y)\in \R \times \R^n$ where
\begin{align*}
f^{\sigma}_{L_\eps, \theta}(x)= W_{L_\eps} x \quad \forall x
\in \R^{N^\eps_{L-1}}
\end{align*}
as well as
\begin{align*}
f^{\sigma}_{i, \theta}(x)= \sigma\left(W_{i1}x+b_i\right)+ W_{i2}x \quad \forall x
\in \R^{N^\eps_{i-1}},~i=1,\,\dots,\, L_\eps-1
\end{align*}
is called the~\emph{realization} of~$\theta$ with~\emph{activation function}~$\sigma$. Here the application of~$\sigma$ is defined to act componentwise i.e. given an index~$i \in \{1,\dots,L_\eps-1\}$ and~$x \in \R^{N^\eps_i}$ we set
\begin{align*}
\sigma(x)=\left(\sigma(x_1),\,\dots,\, \sigma\left(x_{N^\eps_i}\right)\right)^\top.
\end{align*}
By construction,~$V^\eps_\theta$ satisfies the terminal condition
\begin{align*}
V^\eps_\theta(T,y)=\frac{1}{2} \left|Q_2(y-y_d(T))\right|^2 \quad \forall y\in
\R^n.
\end{align*}
Moreover Assumption~\ref{ass:approxsmoothness} is fulfilled as confirmed by the following result.

\begin{theo}\label{thm:approxbynetworks}
For every~$\eps>0$ there exists architectures~$\operatorname{\clR}_\eps$ and~$\theta_\eps\in
\clR_\eps$ such that~$V^\eps \in \clC ^4(\clR_\eps
\times \R \linebreak\times \R^n)$ and~$V^\eps_{\theta_\eps}$ satisfies~\eqref{eq:approxcapa}.
\end{theo}
\begin{proof}
Let us set $h(t,y)= V^*(t,y)$ for $(t,y) \in I \times \bar B_{2\widehat M}(0)$. Then $h$ is twice continuously differentiable on $I
\times \bar B_{2\widehat M}(0)$ and $h(T,y)= \frac{1}{2} |Q_2(y-y^T_d)|^2$. A consequence of the universal approximation theorem implies that for all $\eps >0 $ there exists $\tilde h_\eps \in \M_{\net}$ such that
\begin{equation}\label{eq:aux10}
\left\|h-\tilde h_\eps \right\|_{C^2\left(I\times \bar B_{2\widehat M}(0)\right)} \le
\frac{\varepsilon}{2},
\end{equation}
where $\M_{\net}= \spanl \{\sigma(\tilde w\cdot x+ \tilde b):
\tilde w \in \R^{n+1},\, \tilde b\in\R \}$, see eg~\cite[Theorem~4.1]{P99}, \cite{hornik91}. Let us observe that $\tilde h_\eps$ can be expressed as a residual network. Indeed, since
\[
\tilde h_\eps = \sum_{i=1}^M \tilde c_i \sigma\left(\tilde w_i \cdot x +\tilde b_i\right)
\]
for some $M\in \mathbb{N}$, $\tilde w_i \in \R^{n+1}, \tilde b_i, \tilde c_i \in \R$, choosing $L_\epsilon=2, W_{11}\in \R^{M\times (n+1)}$ with rows $\{\tilde w_i\}_{i=1}^M$,
\begin{equation*}
b_1=\col\left(\tilde b_1,\,\dots,\,\tilde b_M\right), \,W_2=\left(\tilde c_1,\,\dots,\,\tilde c_M\right), \, W_{12}=0,
\end{equation*}
we have $\tilde h_\eps = f^\sigma_{2,\theta}\circ f^\sigma_{1,\theta}$. Moreover, $\tilde h_\eps\in C^4(I\times \bar B_{2\widehat M})$. Following~\eqref{eq:shift} we define
\begin{align*}
V^\eps_{\theta_\eps}(t,y)= \frac{1}{2} \left|Q_2\left(y-y^T_d\right)\right|^2 + \tilde h_\eps(t,y) -
\tilde h_\eps(T,y)\in C^4\left(I\times \bar B_{2\widehat M}\right).
\end{align*}
and estimate
\begin{align*}
\left\| V^\eps_{\theta_\epsilon}(t,y) - V^*(t,y)\right\|_{C^2} &=
\left\| \tilde h_{\eps}(t,y) - \tilde h_{\eps}(T,y) + V^*(T,y) - V^*(t,y)\right\|_{C^2} \\
& \leq 2 \left\|h-\tilde{h}_\eps\right\|_{C^2} \leq \eps,
\end{align*}
where all norms are taken over $I\times \bar B_{2\widehat M}(0)$. This ends the proof.
\end{proof}

\subsection{Piecewise polynomials} \label{subsec:polynomials}
Fix~$\eps_0>0$, and let $\eps\in (0, \eps_0]$ be arbitrarily fixed. Throughout this subsection we assume~$(\mathbf{A.2})$ and in particular we shall make use of the global Lipschitz continuity of $D^2 V^*$ on $\bar K= \bar I \times \bar B_{2\widehat M}(0)$. Since $\bar K$ is compact and hence totally bounded, there exist $n_\eps \in \mathbb{N}$ and $\{(\bar{t}_i,\bar{y}^i_0)\}_{i=1}^{n_\eps}\in \R^{n+1}$ such that
\begin{align*}
\bar{K} \subset \bigcup^{n_\eps}_{i=1} K_i \quad \text{where} \quad K_i= B_{\eps}\left(\left(t_i,\bar{y}^i_0\right)\right).
\end{align*}
Note that we do not highlight the dependence of $(t_i,\bar{y}^i_0)$ and $K_i$ on $\eps$. For each~$i$ define the parametrized polynomial
\begin{align*}
V^\eps_i(A,b,c,t,y)= \left(t-\bar{t}_i, y-\bar{y}^i_0\right)^\top A\left(t-\bar{t}_i, y-\bar{y}^i_0\right)+ b^\top \left(t-\bar{t}_i, y-\bar{y}^i_0\right)+c
\end{align*}
with
\begin{align*}
(A,b,c,t,y) \in \Symm(n+1) \times \R^{n+1} \times \R \times \R^{n+1},
\end{align*}
where ${Sym}(n)$ denotes the space of real symmetric $ n \times n$ matrices. Note that~$V^\eps_i $ is infinitely many times differentiable in all of its arguments.

For each $\eps\in (0,\eps_0]$ we define a special partition of unity $\{\varphi_i\}_{i=1}^{n_\eps}$ subordinate to~$K_i$ with $\varphi_i \colon \R \times \R^n \to [0,1]$, satisfying~$\clC ^4$ and
\begin{equation}\label{eq:kk10}
\begin{cases}
&\supp \varphi_i = \bar K_i, \qquad~\sum^{n_\eps}_{i=1} \varphi_i(t,y)=1, \, \forall (t,y) \in \bar K, \\[1.7ex]
&\left\|D^{j}\varphi_i\right\|_{C\left(\bar K_i\cap \bar K\right)} =\bar \mu \eps^{-j}, \forall i=1,\dots, n_\eps, \; \text{ and } j\in\{1,2\},\\[1.7ex]
&\text{card } \left\{i: \varphi_i(t,y) \neq 0 \right\} \le \frkm \quad \forall (t,y)\in \bar K, \eps\in (0,\eps_0],
\end{cases}
\end{equation}
with $\bar \mu$ and $\frkm$ positive constants independent of $i, (t,y)\in \bar K, \eps\in (0,\eps_0]$. Finally we define
\begin{align*}
\clR_\eps= \bigtimes^{n_\eps}_{i=1} \left(\Symm(n+1) \times \R^{n+1} \times \R\right),
\end{align*}
and introduce the family of parameterized functions on $\R^{n+1}$ by
\begin{align}\label{eq:taylorfunctionchoice}
V^\eps_\theta (t,y)= \frac{1}{2} \left|Q_2(y-y_d(T))\right|^2+ \sum^{n_\eps}_{i=1} \varphi_i(t,y)\left(V^\eps_i(A_i,b_i,c_i,t,y)-V^\eps_i(A_i,b_i,c_i,T,y) \right)
\end{align}
for~$\theta=(A_1,b_1,c_1,\,\dots,\, A_{n_\eps},b_{n_\eps},c_{n_\eps})\in\clR\eps$. Obviously we have~$V^\eps_{\cdot}(\cdot) \in\clC ^4(\clR_\eps \times \R \times \R^n)$ and
\begin{align*}
V^\eps_\theta (T,y)=\frac{1}{2} \left|Q_2\left(y-y^T_d\right)\right|^2 \quad \forall y \in \R^n.
\end{align*}
Thus the final time condition in the HJB equation is fulfilled. Next we show that $V^\eps_{\theta}$ satisfies the approximation property in Assumption~\ref{ass:approxsmoothness} for the particular choice of
\begin{multline}\label{eq:taylorparachoice}
\theta_\eps = \left(\partial_{yy} V^*\left(\bar{t}_1,\bar{y}^1_0\right), \partial_y V^*\left(\bar{t}_1,\bar{y}^1_0\right),V^*\left(\bar{t}_1,\bar{y}^1_0\right), \right.\\
\dots,\left. \partial_{yy} V^*\left(\bar{t}_{n_\eps},\bar{y}^{n_\eps}_0\right), \partial_y V^*\left(\bar{t}_{n_\eps},\bar{y}^{n_\eps}_0\right),V^*\left(\bar{t}_{n_\eps},\bar{y}^{n_\eps}_0\right)\right),
\end{multline}
i.e. $V^\eps_i$ in~\eqref{eq:taylorfunctionchoice} are chosen with
\begin{equation}\label{eq:kk11}
\left(\bar{A}_i,\bar{b}_i,\bar{c}_i\right)=\left(\partial_{yy} V^*\left(\bar{t}_i,\bar{y}^i_0\right), \partial_y V^*\left(\bar{t}_i,\bar{y}^i_0\right),V^*\left(\bar{t}_i,\bar{y}^i_0\right)\right) \quad i=1,\,\dots,\,n_\eps.
\end{equation}

\begin{theo}\label{thm:Taylor}
Let~$V^\eps$ and~$\theta_\eps$ be chosen according to~\eqref{eq:taylorfunctionchoice} and~\eqref{eq:taylorparachoice}, respectively, and suppose that~$(\mathbf{A.2})$ and~\eqref{eq:kk10} are satisfied. Then Assumption~\ref{ass:approxsmoothness} holds.
\end{theo}
\begin{proof}
We already argued that~$V_\theta^\eps$ has the desired regularity. It remains to prove the required approximation capabilities. For abbreviation set~$V^\eps_i(t,y)=V^\eps_i({\bar A_i},{\bar b_i},{\bar c_i},t,y)$, with $({\bar A_i},{\bar b_i},{\bar c_i})$ as in~\eqref{eq:kk11}.


Since~$V^\eps_i$ is the second order Taylor expansion of~$V$ at~$(\bar{t}_i,\bar{y}^i_0)$ we conclude that
\begin{equation}\label{eq:kk12}
\left\|V^*-V^\eps_i\right\|_{C^{2-j}\left(\bar K_i \cap\bar K\right)} \leq \bar c\eps^{j+1}, \quad \text{for } j\in{0,1,2,}
\end{equation}
for some~$\bar c>0$ depending on the global Lipschitz constant of $V$ on $\bar K$, and independent of~$\eps \in (0,\eps_0]$ and~$i$. Still recall that the sets $K_i$ depend on $\eps$. To estimate $V^*(t,y) -V^\eps_\theta(t,y)$ we recall that $V^*(T,y)=\frac{1}{2} |Q_2(y-y^T_d)|^2$, and express $V^*(t,y)$ as $V^*(t,y)= V^*(T,y)+ V^*(t,y)- V^*(T,y)$. This leads to
\begin{multline*}
V^*(t,y)-V^\eps_\theta(t,y) \\
= \sum_{i\,\in\,\left\{1,\,\dots,\,n_\eps\right\}} \varphi_i(t,y)\left(V^*(t,y)-V^\eps_i(t,y)\right) +\sum_{i\,\in\,\left\{1,\,\dots,\,n_\eps\right\}}\varphi_i(T,y)\left(V^*(T,y)-V^\eps_i(T,y)\right),
\end{multline*}
for $(t,y) \in \bar K$. From~\eqref{eq:kk10} and~\eqref{eq:kk10} we deduce that $\|V(t,y) -V^\eps_\theta(t,y)\|_{C(\bar K)}\le 2 \bar c \eps^3$.

For the gradient with respect to~$y$ we proceed similarly. Fixing~$(t,y)\in \bar K $ we estimate
\begin{align*}
\left|\partial_y V^*(t,y)-\partial_y V^\eps_{\theta_\eps}(t,y)\right| \leq D_1+D_2
\end{align*}
where
\begin{align*}
D_1 &=\sum_{i\,\in\,\left\{1,\,\dots,\,n_\eps\right\}} \left\lbrack \varphi_i(t,y) \left|\partial_y V^*(t,y)-\partial_y V^\eps_i(t,y)\right|+ \left| V^\eps_i(t,y)-V^*(t,y)\right| \left|\partial_y \varphi_i(t,y)\right| \right \rbrack \\
D_2 &=\sum_{i\,\in\,\left\{1,\,\dots,\,n_\eps\right\}} \left \lbrack \varphi_i(T,y) \left|\partial_y V^*(T,y)-\partial_y V^\eps_i(T,y)\right|+ \left| V^\eps_i(T,y)-V^*(T,y)\right| \left|\partial_y \varphi_i(t,y)\right| \right \rbrack.
\end{align*}
By~\eqref{eq:kk12} with $j=1$ the first terms in $D_1$ and $D_2$ can be estimated by $\bar c \eps ^2$. Using~\eqref{eq:kk10} and~\eqref{eq:kk12} the second terms in $D_1$ and $D_2$ can be bounded by $\frkm \bar \mu \eps ^2$. Combining these estimate we arrive at
\begin{equation*}
\left\|\partial_y V^*(t,y) -\partial_y V^\eps_\theta(t,y)\right\|_{C(\bar K)}\le 2 \eps^2(\bar c + \frkm\, \bar \mu).
\end{equation*}
In an analogous manner one can obtain a bound of the order $O(\eps)$ on the difference of the Hessians of $V$ and $V^\eps_\theta$. This finishes the proof of Theorem~\ref{thm:Taylor}.
\end{proof}

In Appendix~\ref{app1} it is shown how standard mollifiers can be used so that~\eqref{eq:kk10} is satisfied. This requires some extra attention due to the required bounds on the derivatives of $\varphi_i$.


\section{Existence of minimizers to~\eqref{def:approxfeedprop}} \label{sec:existence}
This section is devoted to proving the existence of minimizing triples to~\eqref{def:approxfeedprop}. Throughout this section $c$ will denote a generic constant independent of $\eps>0$ and $y_0\in Y_0$.

\subsection{Existence of admissible points} \label{subsec:exstofmin}
Recall from Assumption~\ref{ass:feedbacklaw} \tcb and Remark~\ref{rem3} that the optimal ensemble state~$\bfy^* \in \clC (Y_0;W_T)$ satisfies~$\|{\bfy^*}\|_{{\mathcal C}}\leq M_{Y_0}$. Accordingly we define the set of admissible states and admissible controls as
\begin{align*}
\bfy_{ad} = \left\{\, \bfy \in \clC (Y_0;W_T)\;\middle|\;\Cbochnorm{\bfy} \leq 2 M_{Y_0}\, \right\}, \quad \bfU_{ad} \coloneqq L^2\left(Y_0;L^2\left(I;\R^m\right)\right).
\end{align*}
We also recall the definition $\clY _{ad}$ in~\eqref{eq:kk13}.

To prove the existence of minimizers to~\eqref{def:approxfeedprop} we first argue that the admissible set
\begin{align}\label{def:admissibleset}
\clN^{\eps}_{ad}= \left\{(\bfy,\bfp,\theta)\!\in\! \bfy_{ad}\!\times\!\clC (Y_0;W_T) \!\times\! \clR_{\eps}\middle|(\bfy,\bfp,\theta)~\text{satisfies}\!~\eqref{eq:statepropapprox}\! -\!\eqref{eq:constraintprop}, \F^\eps_{\theta}(\bfy)\!\in\! \bfU_{ad}\right\}
\end{align}
is nonempty for~$\eps$ small enough. For this purpose consider the family~$\theta_\eps \in \clR_{\eps}$,~$0< \eps \leq \eps_0$, from Assumption~\ref{ass:approxsmoothness} as well as the associated closed loop system of state and adjoint equations
\begin{align}\label{eq:neuralnetstate}
\dot{y}_\eps&= \bff(y_\epsilon)+ \bfg(y_\eps)
\F^\eps_{\theta_\epsilon} (y_\epsilon), \\
-\dot{p}_\eps&= D\bff(y_\eps)^\top p_\eps + \left\lbrack D\bfg(y_\eps)^\top\F^\eps_{\theta_\eps}(y_\eps)\right\rbrack p_\eps+\bfQ_1^\top \bfQ_1(y_\eps-y_d), \label{eq:neuraladjoint}
\end{align}
subject to the following initial and terminal conditions
\begin{align*}
y_\epsilon(0)=y_0,~p_\eps(T)= Q^\top_2 Q_2\left(y_\eps(T)-y_d^T\right),
\end{align*}
for every~$y_0 \in Y_0$. We first prove the following approximation result.

\begin{theo}\label{thm:existenceneuralnetwork2}
Let Assumptions~\ref{ass:feedbacklaw} and~\ref{ass:approxsmoothness} hold. There exists a constant $c$ such that for all~$\eps>0$ small enough and for all $y_0\in Y_0$ the system~\eqref{eq:neuralnetstate} and~\eqref{eq:neuraladjoint} admits unique solutions~$y_\eps=\bfy_\eps (y_0) \in \clY _{ad}$ and~$p_\eps=\bfp_\eps (y_0) \in W_T$. Furthermore~$\bfy_\eps \in \clC ^1 (Y_0;W_T)$,~$\bfp_\eps \in \clC  (Y_0;W_T)$, and $\F^*(\bfy^*) \in \clC  (Y_0;L^2(I;\R^m))$ hold and
\begin{align*}
\left\|\bfy_\eps-\bfy^*\right\|_{\clC^1(Y_0;W_T)}+
\left\|\bfp_\eps-\bfp^*\right\|_{\clC (Y_0;W_T)}+ \left\|\F^\eps_{\theta_\eps}(\bfy_\eps)-\F^*(\bfy^*)\right\|_{\clC\left(Y_0;L^2(I;\R^m)\right)} \leq c \eps.
\end{align*}
In particular,~$(\bfy_\eps,\bfp_\eps,\theta_\eps)\in \clN^\eps_{ad}$ for all~$\eps>0$ small enough.
\end{theo}

In order to prove this we require several auxiliary results.

\begin{lemm}\label{lem:lipschitzoffeed}
There exists a constant $c$ such that for all $\eps$ small enough there holds
\begin{align*}
\left\|\left(\F^*(y_1)-\F^\eps_{\theta_\eps}(y_1)\right)-\left(\F^*(y_2)-\F^\eps_{\theta_\eps}(y_2)\right)\right\|_{L^2\left(I;\R^m\right)}\leq c\eps \wnorm{y_1-y_2}, \quad \forall y_1,y_2 \in \clY _{ad}.
\end{align*}
\end{lemm}
\begin{proof}
According to the definition of~$\F^*$ and~$\F^\eps_{\theta_\eps}$ we split
\begin{align*}
\left\|\left(\F^*(y_1)-\F^\eps_{\theta_\eps}(y_1)\right)-\left(\F^*(y_2)-\F^\eps_{\theta_\eps}(y_2)\right)\right\|_{L^2\left(I;\R^m\right)} \leq D_1+D_2
\end{align*}
with
\begin{align*}
D_1&= 1/\beta \, \left\|\bfg(y_1)^\top\right\|_{\clB \left(L^2\left(I;\R^n\right),L^2\left(I;\R^m\right)\right)} \left\|\partial_y\left(\left(\clV^*(y_1)\!-\!\clV^\eps_{\theta_\eps}(y_1)\right)-\left(\clV^*(y_2)\!-\!\clV^\eps_{\theta_\eps}(y_2)\right)\right)\right\|_{L^2\left(I;\R^n\right)} \\
D_2 &= 1/\beta \, \left\|\bfg(y_1)^\top-\bfg(y_2)^\top\right\|_{\clB \left(L^2(I;\R^n),L^2(I;\R^m)\right)} \left\|\partial_y\left(\clV^*(y_2)-\clV^\eps_{\theta_\eps}(y_2)\right)\right\|_{L^2\left(I;\R^n\right)}.
\end{align*}
Applying the integral mean value theorem yields
\begin{multline*}
\left\|\partial_y\left(\left(\clV^*(y_1)-\clV^\eps_{\theta_\eps}(y_1)\right)-\left(\clV^*(y_2)-\clV^\eps_{\theta_\eps}(y_2)\right)\right)\right\|_{L^2(I;\R^n)}\\
\leq\sup_{s\,\in\,[0,1]}\left\|\partial_{yy}\left(\clV^*(y_1+sh)-\clV^\eps_{\theta_\eps}(y_1+sh)\right)\right\|_{\clB\left(W_T,L^2(I;\R^n)\right)}\wnorm{h}
\end{multline*}
with~$h=y_2-y_1\in W_T$. Note that~$y_1+sh \in \clY _{ad}$ for all~$s \in [0,1]$. Thus we can use Assumption~\ref{ass:approxsmoothness} for every~$s \in [0,1]$ and~$\delta y \in W_\infty$ and estimate
\begin{multline*}
\left\|\partial_{yy}\left(\clV^*(y_1+sh)-\clV^\eps_{\theta_\eps}(y_1+sh)\right)\delta y\right\|_{L^2\left(I;\R^n\right)}\\
\begin{aligned}
&\leq \sqrt{\int_0^T\left|\partial_{yy}\left(V^*\left(t,y_1(t)+sh(t)\right)-V^\eps_{\theta_\eps}\left(t,y_1(t)+sh(t)\right)\right)\right|^2_{\R^{n
\times n}}|\delta y(t)|^2\de  t}
\\
&\leq c \eps \|\delta y\|_{L^2(I;\R^n)} \leq \eps c \wnorm{\delta y}.
\end{aligned}
\end{multline*}
Similarly we obtain
\begin{align*}
\left\|\partial_{y}\left(\clV^*(y_2)-\clV^\eps_{\theta_\eps}(y_2)\right)\right\|_{L^2(I;\R^n)}= \sqrt{\int^T_0\left|\partial_{y}\left({V}^*(t,y_2(t))-{V}^\eps_{\theta_\eps}(t,y_2(t))\right)\right|^2~\de t} \leq \sqrt{T}c \eps.
\end{align*}
Last recall that~$\bfg$ is Lipschitz continuous and uniformly bounded on~$\clY _{ad}$. Combining these facts yields the desired statement.
\end{proof}

With the same arguments the following a priori estimate can be obtained. For the sake of brevity its proof is omitted.

\begin{coro}\label{lem:calmness}
There exists a constant $c$ such that for all $\eps$ small enough there holds
\begin{align*}
\left\|\F^*(y)-\F^\eps_{\theta_\eps}(y)\right\|_{L^2\left(I;\R^m\right)}
\leq c \eps \wnorm{y}, \quad \forall y \in \clY _{ad}.
\end{align*}
\end{coro}

Next we establish existence of a unique solution to~\eqref{eq:neuralnetstate} as well as a first approximation result.

\begin{prop}\label{thm:existenceneuralnetwork}
Let Assumptions~\ref{ass:feedbacklaw} and~\ref{ass:approxsmoothness} hold. Then for all~$\eps>0$ small enough there is a unique~$\bfy_\eps \in \clC ^1(Y_0;W_T)$ such that~$y_\eps\coloneqq\bfy_\eps(y_0) \in \clY _{ad}$ satisfies~\eqref{eq:neuralnetstate} for all $y_0 \in Y_0$. Moreover there exists a constant $c$ independent of $\eps$ such that
\begin{align*}
\left\|\bfy^*-\bfy_\eps\right\|_{\clC (Y_0;W_T)}+ \left\|\F^*(\bfy^*)-\F^{\eps}_{\theta_\eps}(\bfy_\eps)\right\|_{\clC \left(Y_0; L^2(I;\R^m)\right)} \leq c\eps.
\end{align*}
In particular we have $\|\bfy_\eps\|_{\clC (Y_0;W_T)} \leq 2 M_{Y_0}$ for all sufficiently small $\eps$.
\end{prop}

\begin{proof}
The proof is based on a fixed-point argument. Let~$y_0 \in Y_0$ be arbitrary but fixed. Define the set
\begin{align*}
\M= \left\{ y \in W_T \;|\;\wnorm{y} \leq \frac{3}{2} M_{Y_0}\,\right\} \subset \clY _{ad}.
\end{align*}
On~$\M$ we consider the mapping~$\clZ \colon \M
\to W_T$, where~$z=\clZ(y)\in \clY _{ad}$ is the unique solution of
\begin{align}\label{eq:auxeqfixpoint}
\dot{z}= \bff(z)+ \bfg(z)\F^*(z)+
\bfg(y)\F^\eps_{\theta_\epsilon}(y)-\bfg(y) \F^*(y),
\quad z(0)=y_0.
\end{align}
It is well-defined since the perturbation function~$v=\bfg(y)\F^\eps_{\theta_\epsilon}(y)-
\bfg(y) \F^*(y)\in L^2(I;\R^n)$ satisfies
\begin{align*}
\|v\|_{L^2}\leq \|\bfg(y)\|_{\clB \left(L^2(I;\R^m),L^2(I;\R^n)\right)}
\left\|\F^*(y)-\F^\eps_{\theta_\epsilon}(y)\right\|_{L^2} \leq \frac{3}{2} c \eps M_{Y_0} \|\bfg(y)\|_{\clB \left(L^2(I;\R^m),L^2(I;\R^n)\right)}
\end{align*}
where we use Corollary~\ref{lem:calmness} and the definition of~$\M$. Hence $\|v\|_{L^2}\leq c \eps$. Here and below $c$ denotes a generic constant which is independent of $y_0\in Y_0$ and all $\epsilon>0$ sufficiently small. We may invoke Proposition~\ref{thm:existspert} and Corollary~\ref{coroll:locallipschitzofstate} from the Appendix, to assert the existence of a unique solution~$z\in \clY _{ad}$ to~\eqref{eq:auxeqfixpoint} with
\begin{align*}
\wnorm{z}& \leq M_{Y_0}+c\|v\|_{L^2} \leq \frac{3}{2} M_{Y_0},\quad \forall y_0 \in Y_0,
\end{align*}
if~$\eps>0$ is chosen small enough. From this we particularly conclude~$\clZ(\M)\subset \M$ for all~$y_0 \in Y_0$ and~$\eps >0$ small. It remains to prove that~$\clZ$ is a contraction. To this end let~$y_1,~y_2 \in
\M$ be given. Applying Corollary~\ref{coroll:locallipschitzofstate} yields the first inequality in
\begin{align*}
\left\|\clZ(y_1)-\clZ(y_2)\right\|_{W_T} &\leq c
\left\|\F^*(y_1)-\F^\eps_{\theta_\epsilon}(y_1)-\F^*(y_2)+\F^\eps_{\theta_\epsilon}(y_2)\right\|_{L^2} \leq c \eps \wnorm{y_1-y_2}
\end{align*}
with a constant~$c>0$ independent of~$y_1, y_2 \in \M$ as well as of~$y_0 \in Y_0$, and $\epsilon$ sufficiently small. The last inequality follows from Lemma~\ref{lem:lipschitzoffeed}. Choosing~$\eps>0$ small enough we conclude that~$\clZ$ admits a unique fixed point~$y_\eps=\clZ(y_\eps) \in W_T$ on~$\M$. Clearly, the function~$ \bfy_\eps(y_0):=y_\eps$ satisfies~\eqref{eq:neuralnetstate},~$y_\eps \in {\M} \subset\clY _{ad}$ as well as
\begin{align*}
\left\|\bfy_\eps(y_0)-\bfy^*(y_0)\right\|_{W_T} =\left\|\clZ(\bfy_\eps(y_0))-\clZ(0)\right\|_{W_T}\leq c\eps \wnorm{y_\eps}
\leq c \eps \frac{3}{2} M_{Y_0},
\end{align*}
and by Corollary~\ref{lem:calmness}
\begin{multline*}
\left\|\F^*\left(\bfy^*(y_0)\right)-\F^\eps_{\theta_\eps}(\bfy_\eps(y_0))\right\|_{L^2} \\
\begin{aligned}
&\leq \|\F^*\left(\bfy^*(y_0)\right)-\F^*\left(\bfy_\eps(y_0)\right)\|_{L^2} + \left\|\F^*\left(\bfy_\eps(y_0)\right)-\F^\eps_{\theta_\eps}(\bfy_\eps(y_0))\right\|_{L^2} \\
& \leq c \left\|\bfy^*(y_0)-\bfy_\eps(y_0)\right\|_{W_T}+c \eps \wnorm{\bfy_\eps(y_0)} \leq c \eps.
\end{aligned}
\end{multline*}
Finally according to Proposition~\ref{thm:existspert} the solution $\bfy_\eps(y_0) $ is unique and the mapping~$\bfy_\eps $ is at least of class~$\clC ^1$.
\end{proof}

Next we estimate the~$W^{1,2}$ difference between~$\bfy_\eps$ and~$\bfy^*$.

\begin{prop}\label{thm:aprioriW12}
The mapping~$\bfy_\eps \in \clC ^1(Y_0;W_T)$ from Theorem~\ref{thm:existenceneuralnetwork} satisfies
\begin{align*}
\left\|\bfy_\eps-\bfy^*\right\|_{\clC ^1(Y_0;W_T)} \leq c\eps
\end{align*}
for~$c>0$ independently of~$\eps$ small enough.
\end{prop}
\begin{proof}
By the previous proposition the estimate is already known for ${\clC ^1(Y_0;W_T)}$ replaced by ${\clC (Y_0;W_T)}$. Now fix~$y_0 \in Y_0$ and~$i \in \{1, \dots,n\}$. By the inverse mapping theorem the partial derivatives of~$\bfy^*$ and~$\bfy_\eps$ at~$y_0$ are given by~$\partial_i \bfy^*(y_0)= T_*(y_0)^{-1}(0,e_i)$,~$\partial_i \bfy_\eps(y_0)= T_\eps (y_0)^{-1}(0,e_i)$. Here,~$e_i$ denotes the i$^{\rm th}$ canonical basis vector in~$\R^n$ and
\begin{align*}
T_*(y_0)^{-1}, T_\eps(y_0)^{-1} \colon L^2\left(I;\R^n\right) \times \R^n \to W_T
\end{align*}
denote the linear continuous inverses of
\begin{align*}
T_*(y_0)\delta y\!=\!\left(
\begin{array}{c}
\dot{\delta y}-D\bff\left(\bfy^*(y_0)\right)\delta y-\left\lbrack D\bfg\left(\bfy^*(y_0)\right)\delta y\right\rbrack\F^*\left(\bfy^*(y_0)\right)-\bfg\left(\bfy^*(y_0)\right)D\F^*\left(\bfy^*(y_0)\right)\delta y\\
\delta y(0)\\
\end{array}
\right)
\end{align*}
and
\begin{align*}
T_\eps(y_0)\delta y\!=\!\left(
\begin{array}{c}
\dot{\delta y}-D\bff\left(\bfy_\eps(y_0)\right)\delta y-\left\lbrack D\bfg\left(\bfy_\eps(y_0)\right)\delta y\right\rbrack\F^\eps_{\theta_\eps}(\bfy_\eps(y_0))-\bfg(\bfy_\eps(y_0))D\F^\eps_{\theta_\eps}(\bfy_\eps(y_0))\delta y\\
\delta y(0)\\
\end{array}
\right).
\end{align*}
Using Gronwall's inequality, we readily verify that
\begin{align}\label{eq:uniformT}
\max \left\{ \left\|T_\eps (y_0)^{-1}(\delta v, \delta y_0)\right\|_{W_T},\left\|T_*(y_0)^{-1}(\delta v, \delta y_0)\right\|_{W_T} \right\} \leq C\left(\|\delta v\|_{L^2(I;\R^n)}+|\delta y_0|_{\R^n}\right)
\end{align}
for all~$ \delta v \in L^2(I,\R^n),~\delta y_0 \in \R^n$,~$y_0 \in Y_0$ and some~$C>0$ independent of~$y_0, \delta v, \delta y_0 $. Now we recall that~$\bfy_\eps(y_0), \bfy^*(y_0)\in \clY _{ad}$ and that~$D \bff, D \bfg,\bfg$ are Lipschitz continuous, and thus in particular bounded, on~$\clY _{ad}$, see Assumption $(\mathbf{A.1})$. Together with~boundedness of $\{ \|\F^*(\bfy(y_0))\|_{L^2}: y_0 \in Y_0\}$, Corollary~\ref{lem:calmness} and Theorem~\ref{thm:existenceneuralnetwork} we conclude
\begin{align}\label{eq:kk14}
\left\|\left(T_*(y_0)-T_\eps(y_0)\right)\delta y\right\|_{L^2\left(I;\R^n\right)\times \R^n} \leq c \eps \wnorm{\delta y} \quad \forall \delta y \in W_T
\end{align}
for some~$c>0$ again independent of~$y_0 \in Y_0$. Recalling that $B^{-1} - A^{-1}=A^{-1}(A-B)B^{-1}$ for invertible bounded linear operators $A$ and $B$, we obtain
\begin{align*}
\left\|\partial_i\bfy_\eps(y_0)-\partial_i \bfy^*(y_0)\right\|_{W_T}&= \left\|T_\eps (y_0)^{-1}(0,e_i)-T_* (y_0)^{-1}(0,e_i)\right\|_{W_T} \\
& \leq C^2 \sup_{\|\delta y\|_{W_T}\,\leq\,1} \left\|\left(T_*(y_0)-T_\eps(y_0)\right)\delta y\right\|_{L^2\left(I;\R^n\right)\times \R^n} \leq c \eps,
\end{align*}
where~$C>0$ is the constant from~\eqref{eq:uniformT}. Since all involved constants are independent of~$y_0 \in Y_0$ we obtain the desired estimate $
\left\|\partial_i \bfy_\eps-\partial_i \bfy^*\right\|_{\clC } \leq c \eps$.
\end{proof}

Next we address the solvability of the adjoint equation~\eqref{eq:adjointpropapprox}.

\begin{prop}\label{prop:solvofadjoint}
There exists a constant $c$ such that for all $\eps$ small enough there exists~$\bfp_\eps \in \clC (Y_0;W_T)$ such that~$p_\eps \coloneqq \bfp_\eps(y_0)\in W_T$ satisfies~\eqref{eq:neuraladjoint} for all~$y_0 \in Y_0$ and
\begin{align*}
\left\|\bfp_\eps-\bfp^*\right\|_{\clC} \leq c \eps.
\end{align*}
\end{prop}

\begin{proof}
Given~$y\in \clY _{ad}$ consider the linear ordinary differential equation
\begin{align*}
-\dot{p}= D\bff(y)p+\left\lbrack D \bfg(y)^\top \F^\eps_{\theta_\eps}(y)\right\rbrack p+\bfQ_1^\top \bfQ_1(y-y_d),~p(T)= Q^\top_2 Q_2\left(y(T)-y^T_d\right).
\end{align*}
It admits a unique solution~$p=P(y)\in W_T$ which is bounded independently of~$y\in \clY _{ad}$. Moreover the mapping~$P \colon W_T \to W_T$ is continuous on~$\clY _{ad}$ in virtue of the Gronwall lemma and Assumption~\ref{ass:feedbacklaw}. The existence of a mapping~$\bfp_\eps$ which satisfies~\eqref{eq:neuraladjoint} then follows by setting~$\bfp_\eps=P \circ \bfy_\eps$.

It remains to prove the estimate for the difference between~$\bfp_\eps$ satisfying~\eqref{eq:neuraladjoint} and~$\bfp^*$ satisfying~\eqref{eq:adjointprop}. For this purpose we can use the same technique as in the proof of Proposition~\ref{thm:aprioriW12} and therefore we only give the main estimates. Recall that~$D \bff(\cdot)^\top,D \bfg(\cdot)^\top$ are Lipschitz continuous on~$\clY _{ad}$. The the most involved term in the estimate analogous to~\eqref{eq:kk14} is
\begin{multline*}
\left\|\left\lbrack D \bfg\left(\bfy_\eps(y_0)\right)^\top \F^\eps_{\theta_\eps}\left(\bfy_\eps(y_0)\right)-D \bfg\left(\bfy^*(y_0)\right)^\top \F^*\left(\bfy^*(y_0)\right)\right\rbrack \delta p\right\|_{L^2}\\
\leq c\left(\left\|\bfy_\eps(y_0)-\bfy^*(y_0)\right\|_{W_T}+\left\|\F^\eps_{\theta_\eps}(\bfy_\eps(y_0))-\F^*\left(\bfy^*(y_0)\right)\right\|_{L^2}\right) \wnorm{\delta p}
\end{multline*}
with~$c>0$ independent of~$\eps>0$ and $\delta p \in W_T$. Now a perturbation argument as in the proof of Propostion~\ref{thm:aprioriW12} provides us with
\begin{multline*}
\left\|\bfp_\eps(y_0)-\bfp^*(y_0)\right\|_{W_T}\\
\begin{aligned}
& \leq c \left(\left\|\bfy_\eps(y_0)-\bfy^*(y_0)\right\|_{W_T}+\left\|\F^\eps_{\theta_\eps}\left(\bfy_\eps(y_0)\right)-\F^*\left(\bfy^*(y_0)\right)\right\|_{L^2}+\left|\bfy_\eps(y_0)(T)-\bfy^*(y_0)(T)\right|\right) \\
& \leq c\left(\left\|\bfy_\eps(y_0)-\bfy^*(y_0)\right\|_{W_T}+\left\|\F^\eps_{\theta_\eps}\left(\bfy_\eps(y_0)\right)-\F^*\left(\bfy^*(y_0)\right)\right\|_{L^2}\right) \leq c \eps
\end{aligned}
\end{multline*}
where~$W_T \hookrightarrow \clC (\bar{I};\R^n)$ is used in the second inequality, and Proposition~\ref{thm:aprioriW12} and Corollary~\ref{lem:calmness} are utilized in the final one. Since all involved constants are again independent of~$y_0 \in Y_0$, this finishes the proof.
\end{proof}

Summarizing all previous observations we arrive at the proof of Theorem~\ref{thm:existenceneuralnetwork2}.


\begin{proof}[Proof of Theorem~\ref{thm:existenceneuralnetwork2}]
This follows directly by combining Proposition~\ref{thm:existenceneuralnetwork}, Proposition~\ref{thm:aprioriW12}, and Proposition~\ref{prop:solvofadjoint}.
\end{proof}

\subsection{Closedness of~\texorpdfstring{$\clN^\eps_{ad}$}{N epsad}} \label{subsec:closedandexist}
As a last prerequisite for proving existence to~\eqref{def:approxfeedprop} we argue that the admissible set~$\clN^\eps_{ad}$ is closed. The existence of at least one minimizing triple to~\eqref{def:approxfeedprop} then follows by variational arguments. From here on we always assume that~$\clN^\eps_{ad}$ from~\eqref{def:admissibleset} is nonempty, i.e. that $\eps$ is sufficiently small.



\begin{prop}\label{prop:closedofNad}
Let~$(\bfy_k,\bfp_k,\theta_k)_{k\,\in\,\N} \subset \clN^\eps_{ad} $ be a sequence with weak limit~$(\bfy,\bfp,\theta)$ in~$L^2(Y_0;W_T)^2 \times \clR_\eps$. Then $(\bfy,\bfp,\theta) \in \clN^\eps_{ad} $ and we have
\begin{align*}
(\bfy,\bfp) \in \clC (Y_0;W_T)^2,~\lim_{k\,\rightarrow\,\infty} \bfy_k(y_0)=\bfy(y_0) \text{ and}~\lim_{k\,\rightarrow\,\infty} \bfp_k(y_0)=\bfp(y_0) \text{ in } W_T, \quad \forall~y_0 \in Y_0.
\end{align*}
\end{prop}

The proof builds upon the following two lemmas.

\begin{lemm}\label{lem:weakclosedstate}
Let the sequence $(\bfy_k,\bfp_k,\theta_k)_{k\,\in\,\N} \subset \clN^\eps_{ad} $ satisfy the prerequisites of Proposition~\ref{prop:closedofNad}. Then~$\bfy \in \bfy_{ad}$, $\bfy_k(y_0)\rightarrow \bfy(y_0)\text{ in}~W_T,$ $\F^\eps_{\theta_k}(\bfy_k(y_0))\rightarrow \F^\eps_{\theta}({\bfy}(y_0))~\text{in}~L^\infty (I;\R^m)$, and
\begin{align}\label{eq:aux8}
\dot{\bfy}(y_0)= \bff\left({\bfy}(y_0)\right)+ \bfg\left({\bfy}(y_0)\right)\F^\eps_{\theta}({\bfy}(y_0)),~{\bfy}(y_0)(0)=y_0,
\end{align}
for all~$y_0 \in Y_0$.
\end{lemm}


\begin{proof}
By assumption we have~$\bfy_k \in \bfy_{ad}$, and hence $\|\bfy_k(y_0)\|_{W_T} \le 2M_{Y_0}$ for all $k \in \N$ and $y_0 \in Y_0$, and $\bfy_k \in \clC ^1(Y_0;W_T)$ for all~$k \in \N$, see Proposition~\ref{thm:existenceneuralnetwork}. Let us fix an arbitrary $y_0 \in Y_0$. and set~$y_k \coloneqq \bfy_k(y_0)$ for abbreviation. Then there exists a subsequence, denoted by the same index, and~$\tilde y\in W_T$ such that~$y_k \rightharpoonup \tilde y $ in $W_T$. Since~$W_T \hookrightarrow_c \clC (\bar{I};\R^n) \hookrightarrow L^p(I;\R^n)$,~$1\leq p \leq +\infty$, we immediately get
\begin{align*}
y_k(0) \rightarrow \tilde y(0)~\text{in}~\R^n,~\bff(y_k) \rightarrow \bff(\tilde y)~\text{in}~L^2(I;R^n),~\bfg(y_k) \rightarrow \bfg(\tilde y)~\text{in}~\clB \left(L^2\left(I;\R^m\right),L^2\left(I;\R^n\right)\right)
\end{align*}
as well as $\F^\eps_{\theta}(y_k) \rightarrow \F^\eps_{\theta}(\tilde y)~\text{in}~L^\infty(I;\R^m). $ Moreover by Assumption~\ref{ass:approxsmoothness} for every~$\delta >0$ there exits~$K_\delta\in \N$ such that
\begin{align}\label{eq:kk15}
\left|\partial_y V^\eps_{\theta_k}(t,y)-\partial_y V^\eps_{\theta}(t,y)\right|\leq \delta \quad \forall (t,y)\in \bar{I} \times \bar{B}_{2\widehat{M}}(0)
\end{align}
for all~$k\geq K_\delta$. Here~$\widehat{M}$ denotes the constant from $(\mathbf{A.2})$. For all such~$k$ we get utilizing~\eqref{eq:kk15} for a constant $c$ independent of $k$
\begin{align*}
\left\|\F^\eps_{\theta_k}(y_k)-\F^\eps_{\theta}(\tilde y)\right\|_{L^\infty} &\leq c\left\|\F^\eps_{\theta_k}(y_k)-\F^\eps_{\theta}(y_k)\right\|_{L^\infty}+\left\|\F^\eps_{\theta}(y_k)-\F^\eps_{\theta}(\tilde y)\right\|_{L^\infty} \\
& \leq c \delta + \left\|\F^\eps_{\theta}(y_k)-\F^\eps_{\theta}(\tilde y)\right\|_{L^\infty}.
\end{align*}
This implies that~$\lim_{k\,\to\,\infty} \F^\eps_{\theta_k}(y_k)=\F^\eps_{\theta}(\tilde y)$ in~$L^\infty(I;\R^m)$. These observations imply
\begin{align*}
\dot{y}_k= \bff(y_k)+\bfg(y_k)\F^\eps_{\theta_k}(y_k) \rightarrow \bff(\tilde y)+\bfg(\tilde y)\F^\eps_{\theta}(\tilde y).
\end{align*}
Together with $y_k \rightharpoonup \tilde y $ in $W_T$ this implies that $\bfy_k(y_0)=y_k \rightarrow \tilde y$ in~$W_T$ and
\begin{equation}\label{eq:closedloopaux14}
\dot{\oldtilde{y}}=\bff(\tilde y)+\bfg(\tilde y)\F^\eps_{\theta}(\tilde y), \tilde y(0)=y_0.
\end{equation}
Since the solution to this equation is unique, every weak accumulation point of~$y_k$ satisfies~\eqref{eq:closedloopaux14} and we have $\bfy_k(y_0)\rightarrow \tilde y~\text{in}~W_T$ for the whole sequence. We repeat this construction for all $y_0\in Y_0$. This defines a function $ \tilde {\bfy}:Y_0 \to W_T$ such that $\bfy_k(y_0)\rightarrow \tilde{\bfy}(y_0)~\text{in}~W_T$ and such that~\eqref{eq:closedloopaux14} is satisfied with $\tilde y= \tilde {\bfy}(y_0)$ for each $y_0\in Y_0$. By Proposition~\ref{thm:existenceneuralnetwork} it is the unique solution to~\eqref{eq:aux8}.

Lebesgue's dominated convergence theorem for Bochner integrals~\cite[p.~45]{DU77} implies that $\bfy_k\rightarrow \tilde {\bfy}$ in $L^1(Y_0;W_T)$, and by boundedness of $\{\|\bfy_k\|_{\clC }\}_{k=1}^\infty$ also in $L^2(Y_0;W_T)$. By assumption ${\bfy}_k$ converges weakly in $L^2(Y_0;W_T)$ to $\bfy$. Thus we have $\bfy=\tilde {\bfy}$. Moreover $\|\bfy\|_{\clC } \le 2M_{Y_0}$ and hence $\bfy \in \bfy_{ad}$.
\end{proof}

Next we consider the behavior of the adjoint states~$\bfp_k$.

\begin{lemm}\label{lem:weakclosedadjoint}
Let $(\bfy_k,\bfp_k,\theta_k)_{k\,\in\,\N} \subset \clN^\eps_{ad} $ be a sequence with weak limit~$(\bfy,\bfp,\theta)$ satisfying the prerequisites of Proposition~\ref{prop:closedofNad}. Then~$\|\bfp_k\|_{\clC }\leq C$ for some~$C>0$ and all~$k \in\N$ large enough, and~$\bfp \in \clC (Y_0;W_T)$. Moreover~$\bfp_k(y_0)\rightarrow \bfp(y_0)~\text{in}~W_T$, and
\begin{equation}\label{eq:adjointlim}
\begin{split}
-\dot{\bfp}(y_0) &= D\bff(\bfy(y_0))^\top\bfp(y_0)+ \left[D\bfg(\bfy(y_0))^\top\F^\eps_{\theta}(\bfy(y_0))\right]{\bfp}(y_0)+ \bfQ_1^\top \bfQ_1(\bfy(y_0)-y_d),\\
\bfp(y_0))(T) &= Q^\top_2 Q_2\left(\bfy(T)(y_0)-y^T_d\right),
\end{split}
\end{equation}
for all~$y_0 \in Y_0$.
\end{lemm}
\begin{proof}
From Lemma~\ref{lem:weakclosedstate} recall that for the sequences~$~y_k\coloneqq \bfy_k(y_0) \in \clY _{ad}$ and~$y \coloneqq \bfy(y_0)$ we have for each $y_0\in Y_0$
\begin{align*}
y_k \rightarrow y~\text{in}~W_T,~\F^\eps_{\theta_k}(y_k) \rightarrow \F^\eps_{\theta}(y)~\text{in}~L^\infty(I;\R^m).
\end{align*}
Further for each~$k\in\N$ and~$y_0 \in Y_0 $, the element $p_k \coloneqq \bfp_k(y_0) \in W_T$ satisfies
\begin{align}\label{eq:adjointaux}
-\dot{p}_k= D\bff(y_k)^\top p_k+\left\lbrack D \bfg(y_k)^\top \F^\eps_{\theta_k}(y_k)\right\rbrack p_k+\bfQ_1^\top \bfQ_1(y_k-y_d),~p_k(T)= Q^\top_2 Q_2\left(y_k(T)-y^T_d\right).
\end{align}
Recall from Assumption~\ref{ass:approxsmoothness} that~$\partial_y V_{\cdot}^\eps$ is uniformly continuous on compact sets. Thus for every~$\delta>0$ there is~$K_\delta \in \N$ such that
\begin{align*}
\left|F^\eps_{\theta_k}(t,x)\right| \leq \left|F^\eps_{\theta_k}(t,x)-F^\eps_{\theta}(t,x)\right|+ \left|F^\eps_{\theta}(t,x)\right| \leq \delta+ \max_{(t,x)\,\in\,I \times \bar{B}_{2 \widehat{M}}(0)} \left|F^\eps_{\theta}(t,x)\right| < \infty
\end{align*}
for all~$(t,x)\in I \times \bar{B}_{2 \widehat{M}}(0)$ and~$k\geq K_\delta$. Consequently we obtain
\begin{align*}
\sup_{k\,\ge\,K_\delta} \max_{(t,x)\,\in\,I \times \bar{B}_{2 \widehat{M}}(0)} \|A_k(t,x)\|_{\R^{n\times n}} < \infty, \text{ where } A_k(t,x)= Df(t,x)^{\top}+ D g(t,x)^\top F^\eps_{\theta_k}(t,x).
\end{align*}

%~\text{ and }~
Applying Proposition~\ref{prop:s} to the time-reversed equation~\eqref{eq:adjointaux} implies that
\begin{align*}
\wnorm{p_k} \leq c \left(\left\|\bfQ_1^\top \bfQ_1(y_k-{y}_d)\right\|_{L^2}+ \left|y_k(T)-y^T_d\right| \right)
\end{align*}
for some~$c>0$ independent of~$y_0 \in Y_0$ and all sufficiently large $k$. Since~$\|\bfy_k\|_{\clC }\leq 2M_{Y_0}$ we finally conclude~$\|\bfp_k\|_{\clC } \leq C$ for some~$C>0$ independent of~$k$ sufficiently large. We are now prepared to pass to the limit in~\eqref{eq:adjointaux}. For this purpose we proceed as in the proof of Lemma~\ref{lem:weakclosedstate} and use
\begin{align*}
D \bff(y_k)+ D \bfg(y_k)^\top \F^\eps_{\theta_k}(y_k) &\rightarrow D \bff(y)+ D \bfg(y)^\top \F^\eps_{\theta}(y)~\text{in}~\clB \left(L^2(Y;\R^n)\right),
\\
\intertext{as well as}
\bfQ_1^\top \bfQ_1(y_k-{y}_d) &\rightarrow \bfQ_1^\top \bfQ_1(y-{y}_d)~\text{in}~L^2\left(I;\R^n\right),
\\
\intertext{and}
Q^\top_2 Q_2(y_k(T)-y^T_d) &\rightarrow Q^\top_2 Q_2\left(y(T)-y^T_d\right)~\text{in}~\R^n
\end{align*}
to show that every weak accumulation point~$\tilde p \in W_T$ of~$p_k$ is in fact a strong accumulation point and satisfies the differential equation in~\eqref{eq:adjointlim}. Since the solution to this equation is unique we get~$p_k \rightarrow \tilde p $ in~$W_T$ for the whole sequence. Finally utilizing~$\|\bfp_k\|_{\clC } \leq C$ and Lebesgue's dominated convergence theorem we conclude~$\tilde p=\bfp(y_0)$ for all $y_0 \in Y_0$.
\end{proof}

\begin{proof}[Proof of Proposition~\ref{prop:closedofNad}]
This is a direct consequence of Lemma~\ref{lem:weakclosedstate} and Lemma~\ref{lem:weakclosedadjoint}.
\end{proof}


\subsection{Existence of minimizers}
Finally we prove the existence of at least one minimizing triplet to~\eqref{def:approxfeedprop}.

\begin{theo}\label{thm:existence}
Let Assumption~\ref{ass:feedbacklaw} and~\ref{ass:approxsmoothness} hold. Then for all~$\eps>0$ small enough, Problem~\eqref{def:approxfeedprop} admits at least one minimizing triplet~$(\bfy^*_\eps, \bfp^*_\eps, \theta^*_\eps) \in \clC (Y_0;W_T)^2 \times \clR_\eps $.
\end{theo}

\begin{proof}
According to Theorem~\ref{thm:existenceneuralnetwork2}, the admissible set~$\clN^\eps_{ad}$ is nonempty for~$\eps>0$ small enough. Fix such a~$\eps>0$ and let~$(\bfy_k,\bfp_k, \theta_k) \in \clN^\eps_{ad}$ denote a minimizing sequence for~$\clJ_\eps$ i.e.
\begin{align*}
\clJ_\eps(\bfy_k,\bfp_k, \theta_k) \rightarrow \inf_{(\bfy,\bfp, \theta)\,\in\,\clN^\eps_{ad}} \clJ_\eps(\bfy,\bfp, \theta).
\end{align*}
Since~$\bfy_k \in \bfy_{ad}$ and $
\frac{\gamma_\eps}{2} \|\theta_k\|^2_{\clR_\eps} \leq \clJ_\eps(\bfy_k,\bfp_k, \theta_k), $ for all $k\in\N$, the sequence~$\{(\bfy_k,\theta_k)\}\in L^2(Y_0;W_T)\linebreak\times \clR_\eps$ is bounded. Thus it admits at least one subsequence, denoted by the same index, with
\begin{align*}
(\bfy_k,\theta_k) \rightharpoonup (\bfy^*_\eps,\theta^*_\eps)~\text{in}~L^2(Y_0;W_T) \times \clR_\eps
\end{align*}
for some~$(\bfy^*_\eps,\theta^*_\eps)$. As in the proof of Lemma~\ref{lem:weakclosedadjoint} we verify that~$\|\bfy_k\|_{\clC } \leq C$ and~$\|\bfp_k\|_{\clC } \leq C$ for some~$C>0$ independent of~$k\in\N$. Consequently, by possibly taking another subsequence we arrive at
\begin{align*}
(\bfy_k,\bfp_k,\theta_k) \rightharpoonup \left(\bfy^*_\eps,\bfp^*_\eps,\theta^*_\eps\right)~\text{in}~L^2(Y_0;W_T)^2 \times \clR_\eps
\end{align*}
for some~$(\bfy^*_\eps,\bfp^*_\eps,\theta^*_\eps)\in \clN^\eps_{ad}$. For the following estimates it will be convenient to recall the augmented functional~$J_\eps$, see~\eqref{eq:aux2}, which arises in the running cost of~\eqref{def:approxfeedprop} in compact form:
\begin{align}\label{eq:aux2a}
J_\eps(y,p,\theta)= J\left(y, \F^\eps_\theta(y)\right)+ \frac{\gamma_1}{2} \left\|\clV(y)-J_\bullet \left(y, \F^\eps_{\theta^*_\eps}(y)\right)\right\|^2_{L^2(I;\R)}+\frac{\gamma_2}{2}\left\|p-\partial_y \clV(y)\right\|^2_{L^2\left(I;\R^n\right)},
\end{align}
where $J_t$ was defined below~\eqref{def:valuefunc}. Now fix an arbitrary~$y_0\in Y_0$ and set
\begin{align*}
y_k \coloneqq \bfy_k(y_0),~p_k \coloneqq \bfp_k(y_0),~{y^*} \coloneqq \bfy^*_\eps (y_0),~p \coloneqq \bfp^*_\eps (y_0).
\end{align*}
From Lemma~\ref{lem:weakclosedstate} and Lemma~\ref{lem:weakclosedadjoint} we get
\begin{align*}
y_k \rightarrow \tilde y,~p_k \rightarrow p~\text{in}~W_T,~\F^\eps_{\theta_k}(y_k) \rightarrow \F^\eps_{\theta^*_\eps}(\tilde y)~\text{in}~L^2\left(I;\R^n\right)
\end{align*}
and, again using the uniform continuity of~$V_{\bullet}^\eps$ and~$\partial_y V_{\bullet}^\eps$, we conclude
\begin{align*}
\clV^\eps_{\theta_k}(y_k) \rightarrow \clV^\eps_{\theta^*_\eps}(\tilde y)~\text{in}~L^2(I),~\partial_y \clV^\eps_{\theta_k}(y_k) \rightarrow \partial_y \clV^\eps_{\theta^*_\eps}(\tilde y)~\text{in}~L^2\left(I;\R^n\right),
\end{align*}
as well as the uniform boundedness of~$\clV^\eps_{\theta_k}(\bfy_k)$ and~$\partial_y \clV^\eps_{\theta_k}(\bfy_k)$ in~$\clC (Y_0;L^2(I))$ and $\clC (Y_0;L^2(I;\R^n))$, respectively. Moreover we readily verify that
\begin{align*}
\left|J_t(y_k, \F^\eps_{\theta_k}(y_k))\!-\!J_t\left(\tilde y, \F^\eps_{\theta^*_\eps}(\tilde y)\right)\right|\! \leq\! c\!\left(\left\|y_k-\tilde y\right\|_{L^2}+\left\|\F^\eps_{\theta_k}(y_k)\!-\!\F^\eps_{\theta^*_\eps}(\tilde y)\right\|_{L^2}+ \left|y_k(T)- \tilde y(T)\right| \right),
\end{align*}
for some~$c>0$ independent of~$y_0 \in Y_0$, $t\in (0,T)$, and~$k\in \N$. Thus we arrive at
\begin{align*}
J_{{\bullet}} \left(y_k, \F^\eps_{\theta_k}(y_k)\right) \rightarrow J_\bullet \left(y, \F^\eps_{\theta^*_\eps}(y)\right)~\text{in}~L^\infty(I).
\end{align*}
Summarizing the previous findings there holds
\begin{align*}
\left\|\clV(y_k)\!-\!J_\bullet \left(y_k, \F^\eps_{\theta_k}(y_k)\right)\right\|^2_{L^2}\!+\!\left\|p_k\!-\!\partial_y \clV(y_k)\right\|^2_{L^2} &\to \left\|\clV(y)\!-\!J_\bcdot \left(y, \F^\eps_{\theta^*_\eps}(y)\right)\right\|^2_{L^2}+\left\|p\!-\!\partial_y \clV(y)\right\|^2_{L^2}
\\
J\left(y_k,\F^\eps_{\theta_k}(y_k)\right) & \rightarrow J\left(y,\F^\eps_{\theta^*_\eps}(y)\right).
\end{align*}
Using these expressions in~$J_\eps$ as given in~\eqref{eq:aux2a}, and the boundedness of~$\|y_k\|_{L^2},|y_k(0)|,$ $\|p_k\|_{L^2}~\|\F^\eps_{\theta_k}(y_k)\|_{L^2},~\|\clV^\eps_{\theta_k}(y_k)\|_{L^2}$ independent of~$k \in \N$ and~$y_0 \in Y_0$ we finally get by using Lebesgue's dominated convergence theorem
\begin{equation*}
\clJ_{\eps}\left(\bfy_k,\bfp_k, \theta_k\right) \rightarrow \clJ_{\eps}\left(\bfy^*_\eps,\bfp^*_\eps, \theta^*_\eps\right)= \inf_{(\bfy,\bfp, \theta)\,\in\,\clN^\eps_{ad}} \clJ_\eps(\bfy,\bfp, \theta).\qedhere
\end{equation*}
\end{proof}


\section{Convergence towards optimal controls} \label{sec:convergence}
In Proposition~\ref{thm:existenceneuralnetwork} and~\ref{prop:solvofadjoint} it was established that the ensemble triple~$(\bfy^*, \F(\bfy^*), \bfp^*)$ can be approximated by ensemble triples~$(\bfy_\eps, \F^\eps_{\theta_\eps}(\bfy_\eps), \bfp_\eps)$ in the order $O(\eps)$. In this section, the convergence of solutions to~\eqref{def:approxfeedprop} as~$\eps \rightarrow 0$ is addressed. We first consider the terms in the definition~$\clJ_\eps$, see~\eqref{eq:aux2}. To obtain the desired asymptotic behavior a smallness condition on the regularisation parameter $\gamma_\epsilon$, in relation to the norm of the parameters $\theta_\eps$ describing the approximation quality, is required.

\begin{theo}\label{thm:convofobj}
Let Assumptions~\ref{ass:feedbacklaw} and~\ref{ass:approxsmoothness} hold the latter with~$\theta_\eps \in \clR_\eps$, and let~$(\bfy^*_\eps, \bfp^*_\eps,\theta^*_\eps)$, denote an optimal triple to~\eqref{def:approxfeedprop} for all $\eps>0$ small enough. If additionally $\gamma_\eps \|\theta_\eps\|^2_{\clR_\eps} = O(\eps)$, then
\begin{align*}
0 \leq \int_{Y_0} \omega(y_0) \left \lbrack J\left(\bfy^*_\eps(y_0)\right),\F^\eps_{\theta^*_\eps}\left(\bfy^*_\eps(y_0)\right)- V^*(0,y_0) \right \rbrack~\de  \clL(y_0) \leq c\, \eps
\end{align*}
holds and, if~$\gamma_1,\gamma_2 >0$, we also have
\begin{align*}
\int_{Y_0}\!\! \omega(y_0) \left(\left\| V^\eps_{\theta^*_\eps}\left(t,\bfy^*_\eps(y_0)\right)\!-\!J_\bullet \left(\bfy^*_\eps(y_0), \F^\eps\left(\bfy^*_\eps(y_0)\right)\right)\right\|^2_{L^2}+ \left\|\partial_y V^\eps_{\theta^*_\eps}\left(t,\bfy^*_\eps(y_0)\right)\!-\!\bfp^*_{\eps}(y_0)\right\|^2_{L^2}\right) \de  \clL(y_0)  \leq c \, \eps
\end{align*}
for some~$c>0$ independent of~$\eps$.
\end{theo}
\begin{proof}
Let~$\bfy_\eps,\bfp_\eps$ denote the ensembles of state and adjoint trajectories associated to~$\theta_\eps$, see Theorem~\ref{thm:existenceneuralnetwork2}, for~$\eps>0$ small enough. Then we have
\begin{multline*}
\left| J\left(\bfy_\eps(y_0)\right),\F^\eps_{\theta_\eps}\left(\bfy_\eps(y_0)\right)- V^*(0,y_0) \right|
\\
\leq C \left(\left\|\bfy_\eps(y_0)-\bfy^*(y_0)\right\|_{W_T}+ \left\|\F^\eps_{\theta_\eps}\left(\bfy_\eps(y_0)\right)-\F^*\left(\bfy^*(y_0)\right)\right\|_{L^2(I;\R^m)} \right) \leq c \eps
\end{multline*}
for some~$C>0$ independent of~$\eps$. Here we have used~$V^*(0,y_0)=J(\bfy^*(y_0),\F^*(\bfy^*(y_0)))$ for all~$y_0 \in Y_0$, the embedding~$W_T \hookrightarrow \clC (\bar{I};\R^n)$ as well as the a priori estimates of Proposition~\ref{thm:existenceneuralnetwork}. Next we utilize~$\bfp^*(y_0)=\partial \clV^*(\bfy^*(y_0))$,~$y_0\in Y_0$, to estimate
\begin{multline*}
\left\|\partial_y \clV^\eps_{\theta_\eps}\left(\bfy_\eps(y_0)\right)-\bfp_{\eps}(y_0)\right\|^2_{L^2(I;\R^n)} \\
\begin{aligned}
&\leq 2\left(\left\|\partial_y \clV^\eps_{\theta_\eps}\left(\bfy_\eps(y_0)\right)-\partial_y \clV^*\left(\bfy^*(y_0)\right)\right\|^2_{L^2(I;\R^n)}+\left\|\bfp_{\eps}(y_0)-\bfp^*(y_0)\right\|^2_{L^2(I;\R^n)}\right) \\
&\leq c\eps^2,
\end{aligned}
\end{multline*}
where the last inequality is deduced from Proposition~\ref{thm:existenceneuralnetwork} and Proposition~\ref{prop:solvofadjoint}. Proceeding analogously and using~$V^*(t,\bfy^*(y_0)(t))=J_t(\bfy^*(y_0),\F^*(\bfy^*(y_0)))$ for all~$y_0 \in Y_0$,~$t\in I$, we obtain
\begin{align*}
\int_{Y_0} \omega(y_0) \int^T_0 \left|V^\eps_{\theta_\eps}\left(t,\bfy_\eps(y_0)(t)\right)-J_t\left(\bfy_\eps(y_0)\right),\F^\eps_{\theta_\eps} \left(\bfy_\eps(y_0)\right)\right|^2~\de t\,\de \clL(y_0)\leq D_1+D_2,
\end{align*}
where, using Assumption~\ref{ass:approxsmoothness} and again Proposition~\ref{thm:existenceneuralnetwork}
\begin{align*}
D_1 &\coloneqq \int_{Y_0} \omega(y_0) \left\|\clV^\eps_{\theta_\eps}\left(\bfy_\eps(y_0)\right)-\clV^*\left(\bfy^*(y_0)\right)\right\|^2_{L^2(I)}~\de  \clL(y_0) \leq c \eps^2, \\
D_2 & \coloneqq \int_{Y_0}\omega(y_0)\int^T_0 \left|J_t\left(\bfy^*(y_0)\right),\F^*\left(\bfy^*(y_0)\right)-J_t\left(\bfy_\eps(y_0)\right),\F^\eps_{\theta_\eps}\left(\bfy_\eps(y_0)\right)\right|^2~\de t\de \clL(y_0) \leq c \eps^2.
\end{align*}
Combining the previous estimates with the optimality of~$(\bfy^*,\bfp^*,\theta^*_\eps)$, and the assumption on the asymptotic behavior of $\gamma_\epsilon$ we deduce that
\begin{align*}
0 &\leq \int_{Y_0} \omega(y_0)\left \lbrack J_\eps\left(\bfy^*_\eps(y_0),\bfp^*_\eps(y_0),\theta^*_\eps\right)- V^*(0,y_0) \right \rbrack~\de  \clL(y_0) +\frac{\gamma_\eps}{2} \left\|\theta^*_\eps\right\|^2_{\clR_\eps}\\
&\leq \int_{Y_0} \omega(y_0) \left \lbrack J_\eps\left(\bfy_\eps(y_0),\bfp_\eps(y_0),\theta_\eps\right)- V^*(0,y_0) \right \rbrack~\de  \clL(y_0) +\frac{\gamma_\eps}{2} \|\theta_\eps\|^2_{\clR_\eps}
\leq c\, \eps.
\end{align*}
Recalling the definition of $J_\eps$, this yields all claimed estimates and finishes the proof.
\end{proof}

Next the convergence of the ensemble trajectories~$(\bfy^*_\eps, \bfp^*_\eps)$, the feedback controls~$\F^\eps_{\theta^*_\eps}(\bfy^*_\eps)$ as well as the approximate value function~$\clV^\eps_{\theta^*_\eps}$ are analyzed. For this purpose we make use of the additional regularity of ensemble solutions to the closed loop system, see Proposition~\ref{thm:aprioriW12}, and introduce further constraints to~\eqref{def:approxfeedprop}. Without changing the notation we henceforth set
\begin{align}\label{def:Yadstrict}
\mathbf{\hat Y}_{ad} = \left\{\bfy \in \clC ^1(Y_0;W_T)\;\middle|\;\Cbochnorm{\bfy} \leq 2 M_{Y_0},~\|\bfy\|_{W^{1,2}}\leq 2 M_{W^{1,2}}\, \right\},
\end{align}
where~$M_{W^{1,2}}>0$ is a constant with~$\|\bfy^*\|_{W^{1,2}} \leq M_{W^{1,2}} $, the function $y^*$ was introduced in~\eqref{A.3}, and $W^{1,2}= \{\bfy \in L^2(I;W_T): \partial_i \bfy \in L^2(Y_0;W_T), i\in \{1,\,\dots,\,n\} \}$ endowed with the natural norm. Next we note that
\begin{align*}
\frac{\beta}{2} \left\|\F^*\left(\bfy^*(y_0)\right)\right\|^2_{L^2} \leq J\left(\bfy^*(y_0),\F^*\left(\bfy^*(y_0)\right)\right)=V^*(0,y_0)
\end{align*}
for all~$y_0 \in Y_0$. Thus, due to the continuity of the value function~$V^*$, see $(\mathbf{A.2})$, there is~$M_U >0$ with~$\|\F^*(\bfy^*)\|_{L^\infty} \leq M_U$. Correspondingly we set
\begin{align}\label{def:Uadstricter}
\mathbf{\hat U}_{ad} = \left\{\bfu \in L^\infty\left(Y_0;L^2\left(I;\R^m\right)\right)\;|\; \|\bfu\|_{L^\infty}\leq 2 M_U \,\right\}.
\end{align}
We point out that Theorem~\ref{thm:existence} remains valid despite the additional restriction of the set of admissible states and controls. Problem~\eqref{def:approxfeedprop} with $\bfY_{ad}, \bfU_{ad}$ replaced by $\mathbf{\hat Y}_{ad}, \mathbf{\hat U}_{ad}$ will be denoted by $(\mathcal{\hat P}_\epsilon) $.

\begin{prop}\label{prop:existencestrict}
Let Assumption~\ref{ass:feedbacklaw} and~\ref{ass:approxsmoothness} hold. Then for all~$\eps>0$ small enough, Problem $(\mathcal{\hat P}_\epsilon) $ admits at least one minimizing triple.
\end{prop}

\begin{proof}
Let~$(\bfy_\eps,\bfp_\eps, \theta_\eps)$ be defined as in Theorem~\ref{thm:existenceneuralnetwork2}. Then we have~$\bfy_\eps \in \mathbf{\hat Y}_{ad}$, see Proposition~\ref{thm:existenceneuralnetwork} and Proposition~\ref{thm:aprioriW12}, as well as~$\F^\eps_{\theta_\eps}(\bfy_\eps) \in \mathbf{\hat U}_{ad}$, according to Proposition~\ref{thm:existenceneuralnetwork}, for all~$\eps>0$ small enough. Hence the admissible set of~$(\mathcal{\hat P}_\epsilon) $ is not empty. The existence of a minimizing triple then follows by repeating the arguments of the proof of Theorem~\ref{thm:existence} noting that the admissible set
\begin{align*}
\left\{(\bfy,\bfp,\theta)\!\in\! \mathbf{\hat Y}_{ad}\!\times\!\clC (Y_0;W_T) \!\times\! \clR_{\eps}\,\middle|\,(\bfy,\bfp,\theta)~\text{satisfies}\!~\eqref{eq:statepropapprox}\! -\!\eqref{eq:constraintprop}, \F^\eps_{\theta}(\bfy)\!\in\! \mathbf{\hat U}_{ad}\right\}
\end{align*}
is closed w.r.t to the weak topology on~$L^2(Y_0;W_T)^2 \times \clR_\eps$.
\end{proof}

Let us next address the convergence of the optimal ensemble states~$\bfy^*_\eps$, adjoint states~$\bfp_\eps$ and the associated feedback controls~$\F^\eps_{\theta^*_\eps}(\bfy^*_\eps)$ as $\eps$ tends to 0.

\begin{theo}\label{thm:convoftraj}
Let the prerequisites of Theorem~\ref{thm:convofobj} hold, and let $\eps_k >0$ be a strictly decreasing null sequence such that $(\mathcal{\hat P}_{\eps_k})$ admits a minimizing triple $(\bfy^*_k, \bfp^*_k, \theta^*_k)$. Then $(\bfy^*_{k}, \bfp^*_{k}, \F^{\eps_k}_{\theta^*_{\eps_k}}(\bfy^*_{\eps_k}))$ contains at least one accumulation point~$(\bar{\bfy}, \bar{\bfp}, \bar{\bfu}) \in L^\infty(Y_0; W_T)^2 \times L^\infty(Y_0; L^2(I;\R^m))$ w.r.t the strong topology on~$L^2(Y_0;W_T)^2 \times L^2(Y_0;L^2(I;\R^m))$. For each accumulation point and~$\clL$-a.e.~$y_0 \in Y_0$ we have that~$(\bar{y}, \bar{p}, \bar{u}) \coloneqq (\bar{\bfy}(y_0),\bar{\bfp}(y_0),\bar{\bfu}(y_0))$ satisfies
\begin{align*}
(\bar{y}, \bar{u}) \in \min~\eqref{def:openloopproblem}
\end{align*}
as well as
\begin{align*}
\dot{\bar{y}}&=\bff(\bar{y})+\bfg(\bar{y})\bar{u},~\bar{y}(0)=y_0, \\
-\dot{\bar{p}}&= D\bff(\bar{y})^\top \bar{p}+\left\lbrack D \bfg(\bar{y})^\top\bar{u})\right\rbrack \bar{p}+\bfQ_1^\top \bfQ_1(\bar{y}-y_d),~\bar{p}(T)=Q^\top_2 Q_2\left(\bar{y}(T)-y^T_d\right).
\end{align*}
\end{theo}
\begin{proof}
By choice of the admissible sets~$\mathbf{\hat Y}_{ad}$ and~$\mathbf{\hat U}_{ad}$ we have that~$\{(\bfy^*_k, \F^\eps_{\theta^*_k}(\bfy^*_k))\}_{k=1}^{\infty} $ is bounded in~$(W^{1,2}(Y_0;W_T) \cap L^\infty(Y_0;W_T))\times L^\infty(Y_0;L^2(I;\R^m))$. By Gronwall's inequality we can argue that~$\{\bfp^*_k\}_{k=1}^\infty$ is also bounded in~$L^\infty(Y_0;W_T)$. Thus, due to the Banach--Alaoglu theorem, there is a subsequence, denoted by the same index, and~$(\bar{\bfy},\bar{\bfp},\bar{\bfu})\in L^\infty(Y_0;W_T)^2 \times L^\infty(Y_0;L^2(I;\R^m))$ such that
\begin{align*}
\left(\bfy^*_k, \bfp^*_k, \F^{\eps_k}_{\theta^*_k}\left(\bfy^*_k\right)\right) \rightharpoonup^* \left(\bar{\bfy},\bar{\bfp},\bar{\bfu}\right)~\text{in}~L^\infty(Y_0;W_T)^2 \times L^\infty\left(Y_0;L^2(I;\R^m)\right),
\end{align*}
and $\dot{\bfy}^*_k \rightharpoonup \dot{\bar{\bfy}}$ in $L^2(Y_0;L^2(I;\R^n))$. By the compact embedding of $W^{1,2}(Y_0;W_T)$ into $L^2(Y_0;\clC (I;\R^n))$, see~\cite[Theorem~5.3]{AK18} the subsequence can be chosen such that $\bfy^*_k \to \bar{\bfy}$ strongly in $L^2(Y_0;\clC (I;\R^n))$. These properties imply that~$(\bar{y},\bar{u})\coloneqq (\bar{\bfy}(y_0),\bar{\bfu}(y_0))$ satisfies
\begin{equation}\label{eq:aux9}
\dot{\bar{y}}=\bff(\bar{y})+\bfg(\bar{y})\bar{u},~\bar{y}(0)=y_0,
\end{equation}
for~$\clL$-a.e.~$y_0 \in Y_0$. This also implies~$V^*(0,y_0) \leq J(\bar{y},\bar{u}) $ and thus, together with
\begin{align*}
J\left(\bfy^*_k, \F^{\eps_k}_{\theta^*_k}\left(\bfy^*_k\right)\right) \rightarrow V^*(0,\cdot)~\text{in}~L^1(Y_0;W_T), \\
\end{align*}
see Theorem~\ref{thm:convofobj}, we have~$(\bar{y},\bar{u}) \in \argmin~\eqref{def:openloopproblem}$ for~$\clL$-a.e.~$y_0 \in Y_0$. Moreover, again using the strong convergence of $\bfy^*_k$ in $L^2(Y_0;\clC (I;\R^n))$ and recalling the definition of~$J(\cdot,\cdot)$ as
\begin{align*}
J(y,u)= (1/2) \left\|\bfQ_1(y-y_d)\right\|^2_{L^2(I;\R^n)}+(\beta/2) \|u\|^2_{L^2(I;\R^n)}+(1/2) \left|{Q}_2\left(y(T)-y^T_d\right)\right|^2,
\end{align*}
for all~$y \in W_T,~u \in L^2(I;\R^m)$, we also conclude the convergence of the~$L^2(Y_0;L^2(I;\R^m))$ norm of~$\F^{\eps_k}_{\theta^*_k}(\bfy^*_k)$ towards the norm of~$\bar{\bfu}$. Thus~$\F^{\eps_k}_{\theta^*_k}(\bfy^*_k) \rightarrow \bar{\bfu}$ strongly in~$L^2(Y_0;L^2(I;\R^m))$, and $\bfy^*_k \to \bar{\bfy}$ strongly in $L^2(Y_0;W_T)$, by Lebegue's bounded convergence theorem.


It remains to address the strong convergence of~$\bfp_k$. For this purpose we show that the functions~$\lbrack D\bfg(\bfy^*_k)^\top \F^{\eps_k}_{\theta^*_k}(\bfy^*_k(\cdot))\rbrack \bfp^*_k(\cdot) $ converge weakly to~$\lbrack D\bfg(\bar{\bfy}(\cdot))^\top \bar{\bfu}(\cdot)\rbrack \bar{\bfp}(\cdot) $ in~$L^2(Y_0;L^2(I;\R^n))$. Fixing a test function~$\varphi \in L^2(Y_0;L^2(I;\R^n))$ we first note that
\begin{align*}
\lim_{k\,\rightarrow\,\infty}\left(\varphi, \left\lbrack D\bfg(\bar{\bfy}(\cdot))^\top \bar{\bfu}(\cdot))\right\rbrack \left(\bfp^*_k(\cdot)-\bar{\bfp}\right) \right) _{L^2\left(Y_0;L^2(I;\R^n)\right)}=0.
\end{align*}
Second, for~$\clL$-a.e.~$y_0 \in Y_0$ we estimate
\begin{multline*}
\left(\varphi(y_0),\left\lbrack D\bfg\left(\bfy^*_k(y_0)\right)^\top \F^{\eps_k}_{\theta^*_k}\left(\bfy^*_k(y_0)\right)-D\bfg\left(\bar{\bfy}(y_0)\right)^\top \bar{\bfu}(y_0)\right\rbrack \bfp^*_k(y_0)\right)_{L^2\left(I;\R^n\right)} \\
\begin{aligned}
&\leq \! C \|\varphi(y_0)\|_{L^2} \left\|\bfp^*_k(y_0)\right\|_{W_T}\! \left(\left\|\F^{\eps_k}_{\theta^*_k}\left(\bfy^*_k(y_0)\right)\right\|_{L^2}\left\|\bfy^*_k(y_0)\!-\!\bar{\bfy}(y_0)\right\|_{W_T}\!+
\!\left\|\F^{\eps_k}_{\theta^*_k}\left(\bfy^*_k(y_0)\right)\!-\!\bar{\bfu}(y_0)\right\|_{L^2}\!\right) \\
& \leq \! C \|\varphi(y_0)\|_{L^2} \left(\left\|\bfy^*_k(y_0)-\bar{\bfy}(y_0)\right\|_{W_T}+\left\|\F^{\eps_k}_{\theta^*_k}\left(\bfy^*_k(y_0)\right)-\bar{\bfu}(y_0)\right\|_{L^2}\right)
\end{aligned}
\end{multline*}
for some~$C>0$ independent of~$k\in\N$ and~$y_0$. Here we made use of the boundedness of~$\{\bfy^*_k\}_{k=1}^\infty$ and~$\{\bfp^*_k\}_{k=1}^\infty$ in~$L^\infty(Y_0;W_T)$, and of~$\{\F^{\eps_k}_{\theta^*_k}(\bfy^*_k)\}_{k=1}^\infty$ in~$L^\infty(Y_0;L^2(I;\R^m))$. Integrating both sides of the inequality w.r.t to~$\clL$ and utilizing the strong convergence of~$\bfy^*_k$ and~$\F^{\eps_k}_{\theta^*_k}(\bfy^*_k)$ we finally arrive at
\begin{align*}
\lim_{k\,\rightarrow\,\infty}\left(\varphi,\left\lbrack D\bfg\left(\bfy^*_k(\cdot)\right)^\top \F^{\eps_k}_{\theta^*_k}\left(\bfy^*_k(\cdot)\right)-D\bfg(\bar{\bfy}(\cdot))^\top \bar{\bfu}(\cdot)\right\rbrack \bfp^*_k(\cdot)\right)_{L^2\left(Y_0,L^2(I;\R^n)\right)}=0.
\end{align*}
By repeating this argument for the different terms appearing in the adjoint equation we get that~$(\bar{y}, \bar{p}, \bar{u}) \coloneqq (\bar{\bfy}(y_0),\bar{\bfp}(y_0),\bar{\bfu}(y_0))$ satisfies
\begin{align*}
-\dot{\bar{p}}&= D\bff(\bar{y})^\top \bar{p}+\left\lbrack D \bfg(\bar{y})^\top \bar{u}\right\rbrack \bar{p}+\bfQ_1^\top \bfQ_1(\bar{y}-y_d),~\bar{p}(T)= Q^\top_2 Q_2\left(\bar{y}(T)-y^T_d\right)
\end{align*}
for~$\clL$-a.e.~$y_0 \in Y_0$. Applying Gronwall's inequality we deduce
\begin{align*}
\left\|\bfp_k(y_0)-\bar{\bfp}(y_0)\right\|_{W_T} \leq C \left(\left\|\bfy^*_k(y_0)-\bar{\bfy}(y_0)\right\|_{W_T}+\|\F^{\eps_k}_{\theta^*_k}\left(\bfy^*_k(y_0)\right)-\bar{\bfu}(y_0)\|_{L^2} \right)
\end{align*}
for~$\clL$-a.e.~$y_0 \in Y_0$ and~$C>0$ independent of~$y_0$ and~$k$. This yields~$\bfp_k \rightarrow \bar{\bfp}$ strongly in~$L^2(Y_0;W_T)$. Since the weakly convergent subsequence was chosen arbitrarily in the beginning, this finishes the proof.
\end{proof}
\begin{rema}
If~$\bfg(y(t))=B \in \R^{m \times n}$ then the statement of the previous theorem also holds~\text{without} constraints on the control (i.e. for~$\bfU_{ad}=L^2(Y_0;L^2(I;\R^m))$). In this particular case, the uniform boundedness of~$\F^{\eps_k}_{\theta^*_k}(\bfy^*_k)$ in~$L^2(Y_0;L^2(I;\R^m))$ follows from
\begin{align*}
\frac{\beta}{2} \left\|\F^{\eps_k}_{\theta^*_k}\left(\bfy^*_k\right)\right\|^2_{L^2} \leq c \int_{Y_0} \omega(y_0) J\left(\bfy^*_k,\F^{\eps_k}_{\theta^*_k}\left(\bfy^*_k(y_0)\right)\right) \de \clL(y_0) \leq C,
\end{align*}
see Theorem~\ref{thm:convofobj}. Moreover the adjoint equation does no longer depend on the control. Repeating the arguments of the last proof yields the subsequential convergence of~$(\bfy^*_k,\bfp^*_k,\F^{\eps_k}_{\theta^*_k}(\bfy^*_k))$ towards an element
\[
(\bar{\bfy},\bar{\bfp},\bar{\bfu})\in L^\infty(Y_0;W_T)^2 \times L^2\left(Y_0;L^2\left(I;\R^m\right)\right)
\]
such that~$(\bar{y},\bar{p},\bar{u})\coloneqq(\bar{\bfy}(y_0),\bar{\bfp}(y_0),\bar{\bfu}(y_0))$ satisfy the system of state and adjoint equations as well as~$(\bar{y},\bar{u})\in \argmin~\eqref{def:openloopproblem}$ for~$\clL$-a.e.~$y_0 \in Y_0$. Then it only remains to argue the additional regularity~$\bfu\in L^\infty(Y_0;L^2(Y_0;W_T))$. This is, however, a direct consequence of the first order necessary optimality condition~$\bar{\bfu}=(-1/\beta)B^\top\bar{p}$ for~\eqref{def:openloopproblem}, see Proposition~\ref{prop:structure}.
\end{rema}

We point out that the statement of Theorem~\ref{thm:convoftraj} holds independently of the values of the penalty parameters~$\gamma_1,\gamma_2$. If~$\gamma_1,\gamma_2 >0$ then we additionally obtain the following convergence results for the approximate value function~$\clV^\eps_{\theta^*_k}$ and its derivative~$\partial_y \clV^\eps_{\theta^*_k}$ along optimal state trajectories.

\begin{prop}\label{prop18}
Let the prerequisites of Theorem~\ref{thm:convofobj} hold and let~$(\bfy^*_k, \bfp^*_k, \theta^*_k)$ denote a sequence of minimizing triplets as described in Theorem~\ref{thm:convoftraj}. Assume that~$(\bfy^*_{k}, \bfp^*_{k}, \F^{\eps_k}_{\theta^*_k}(\bfy^*_k))$ converges to~$(\bar{\bfy},\bar{\bfp}, \bar{\bfu})$ in
\[
L^2(Y_0;W_T)^2 \times L^2\left(Y_0; L^2\left(I;\R^m\right)\right)
\quad\text{and}\quad\gamma_1,\gamma_2 >0.
\]
Then we also have
\begin{align*}
\clV^{\eps_k}_{\theta^*_k}\left(\bfy^*_k\right) \rightarrow \clV^*(\bar{\bfy})\;\text{in}\;L^2\left(Y_0;L^2(I)\right),\,\partial_y \clV^{\eps_k}_{\theta^*_k}\left(\bfy^*_k\right) \rightarrow \bar{\bfp}~\text{in}~L^2\left(Y_0;L^2\left(I;\R^n\right)\right).
\end{align*}
\end{prop}
\begin{proof}
Due to the convergence of~$\bfy^*_k \to \bar{\bfy}$ in $L^2(Y_0;L^2(I;\R^n))$ and~$\F^{\eps_k}_{\theta^*_k}(\bfy^*_k) \to \bar{\bfu}$ in $L^2(Y_0;L^2(I;\R^m))$, we conclude that
\begin{align*}
J_\bullet \left(\bfy^*_k,\F^{\eps_k}_{\theta^*_k}\left(\bfy^*_k\right)\right) \rightarrow J_\bullet (\bar{\bfy},\bar{\bfu})= \clV^*(\bar{\bfy})~\text{in}~L^2\left(Y_0;L^2(I)\right).
\end{align*}
Together with
\begin{align*}
\lim_{k\,\rightarrow\,\infty} \int_{Y_0} \omega(y_0) \left\| V^{\eps_k}_{\theta^*_k}\left(\bfy^*_k(y_0)\right)-J_\bullet \left(\bfy^*_k(y_0), \F^{\eps_k}\left(\bfy^*_k(y_0)\right)\right)\right\|^2_{L^2}~\de \clL(y_0)=0,
\end{align*}
see Theorem~\ref{thm:convofobj}, we arrive at~$\clV^{\eps_k}_{\theta^*_k}(\bfy^*_k) \rightarrow \clV^*(\bar{\bfy})$ in~$L^2(Y_0;L^2(I))$. The statement on the convergence of~$\partial_y \clV^{\eps_k}_{\theta^*_k}(\bfy^*_k)$ follows similarly from the strong convergence of~$\bfp_k$.
\end{proof}


\section{Learning from a finite training set}\label{sect8}

We turn to analysing a discrete version of~\eqref{def:approxfeedprop}. In this case we can proceed without the state-space constraint $\bfy \in \bfY_{ad}$ provided certain growth bounds on $\bff$ and $\bfg$ are satisfied. The numerical realization of~\eqref{def:approxfeedprop} will always rely on such a discrete approximation. Henceforth we fix a finite ensemble of initial conditions $\{y_0^i: i=1,\,\dots,\,N\} \subset Y_0$. For positive weights $\omega_i$, $i=1,\,\dots,\,N$, and $\eps>0$ we consider
\begin{align}\label{eq:learningprobfinite}
\inf_{y_i, p_i\,\in\,W_T, \theta\,\in\,\clR_\eps} \left
\lbrack\sum^N_{i=1} \omega_i J_\eps(y_i,p_i, \theta)+
\frac{\gamma_\eps}{2} \|\theta\|^2_{\clR_\eps} \right \rbrack
\tag{$\clP^N_\eps$}
\end{align}
subject to
\begin{align*}
\dot{y}_i &=\bff(y_i)+\bfg(y_i)\F^\eps_\theta(y_i),~y_i(0)=y^i_0 \\
-\dot{p}_i&= D\bff(y_i)^\top p_i+ \left\lbrack D\bfg(y_i)^\top\F^\eps_\theta(y_i)\right\rbrack p_i+\bfQ_1^\top
\bfQ_1(y_i-y_d),~p_i(T)= Q^\top_2 Q_2\left(y_i(T)-y_d^T\right).
\end{align*}

Throughout this section, Assumptions~\ref{ass:feedbacklaw} and~\ref{ass:approxsmoothness} are supposed to hold. Further $\eps$ is supposed to be sufficiently small so that the set of admissible solutions for~\eqref{eq:learningprobfinite} is nonempty, compare Theorem~\ref{thm:existenceneuralnetwork2}. It will be convenient to introduce $\bfy=
\col(y_1,\,\dots,\,y_N)$, and $\bfp=
\col(p_1,\,\dots,\,p_N)$, which replace the ensemble states and costates from the previous sections.

\begin{prop}\label{prop:existencefinite}
Let~$\eps>0$ be sufficiently small and let $(\bfy^k, \bfp^k,\theta_k)\in W_T^{2N}\times \clR_\eps$ denote an infimizing sequence for~\eqref{eq:learningprobfinite}. If~$\max_{i} \|y^k_i\|_{L^\infty(I;\R^n)} \leq M_\infty$ for some~${M_\infty}>0$ independent of~$k \in \N$, then Problem~\eqref{eq:learningprobfinite} admits at least one minimizer~$(\bfy^*, \bfp^*, \theta^*)$.
\end{prop}

\begin{proof}
Since by assumption $(\bfy^k, \bfp^k, \theta_k)$ is an infimizing sequence for~\eqref{eq:learningprobfinite} and since~$\beta >0$ we have
\begin{align}\label{eq:dimdepend}
\max_{i} \left\|{\bf{Q_1}} y^k_i\right\|^2_{L^2}+ \max_{i}\left\|\F^\eps_{\theta_k}\left(y^k_i\right)\right\|^2_{L^2} \leq C_N
\end{align}
for some~$C_N >0$ depending on~$N$. Moreover there holds
\begin{align*}
\left\|\dot{y}^k_i\right\|_{L^2} \leq \left\|\bff\left(y^k_i\right)\right\|_{L^2}+
\left\|\bfg(y_i)\F^\eps_{\theta_k}\left(y^k_i\right)\right\|_{L^2} \leq C(\bff,\bfg) M_\infty (1+C_N)
\end{align*}
using the uniform~$L^\infty$ and~$L^2$ boundedness of~$y^k_i$ and~$\F^\eps_\theta(y^k_i)$, respectively. Thus we also have~$\wnorm{y^i_k}\leq \widehat{C}_N $ for all~$k \in \N$, for some~$\widehat{C}_N>0$ which depends on~$N$ but not on~$k$ and~$i$. The proof can now be completed by the same steps as Theorem~\ref{thm:existence}.
\end{proof}

\begin{rema}\label{rem:linftybound}
The~$L^\infty$-boundedness of the minimizing sequence~$y^k_i$ in Proposition~\ref{prop:existencefinite} can be be ensured by additional assumptions on the dynamics of the problem. These include:
\begin{itemize}
\item Add an additional state constraint~$\|y_i\|_{L^\infty} \leq
\widehat{M}$ to~\eqref{eq:learningprobfinite}.
\item Assume that there are~$a_1,a_2,a_3>0$ such that
\begin{align*}
|f(x)| \leq a_1+ a_2 |x|+ a_3 |x|^2,~\|g(x)\| \leq a_1+a_2|x| \quad
\forall x \in \R^n,
\end{align*}
and that $Q_1$ is positive definite. Then by~\eqref{eq:dimdepend} the family $\{y^k_i\}$ is uniformly w.r.t. $i\in\{1,\,\dots,\,n\}$ and $k=1,\,\dots$ bounded in $L^2(I;\R^n)$. Further we can readily verify that
\begin{align*}
\left\|\dot{y}^k_i\right\|_{L^1} &\leq \left\|\bff\left(y^k_i\right)\right\|_{L^1}+
\left\|\bfg\left(y^k_i\right)\right\|_{L^1}
\\ & \leq 2a_1 T+ a_2 \left\|y^k_i\right\|_{L^1}+ a_3 \left\|y^k_i\right\|^2_{L^2}+ a_2
\left\|y^k_i\right\|_{L^2} \left\|\F^\eps_{\theta_k}\left(y^k_i\right)\right\|^2_{L^2} \leq M_N
\end{align*}
for an~$N$-dependent bound~$M_N>0$. Here we made use of the $L^2$-boundedness of~$y^k_i$ and~$\F^\eps_{\theta_k}(y^k_i)$ which follows from~\eqref{eq:dimdepend} in the proof of Proposition~\ref{prop:existencefinite}, and the assumption that $Q_1>0$. Consequently~$y^k_i$ is uniformly bounded in~$W^{1,1}(I; \R^n)$ and thus also in~$L^\infty(I;\R^n)$.
\item Assume that~$f(x)=Ax-h(x)$ where~$A \in \R^{n \times n}$ and~$h$ is monotone i.e.~$(x,h(x))_{\R^n} \geq 0$ for all~$x\in\R^n$. Moreover assume that $Q_1$ is positive definite and that
\begin{align*}
\|g(x)\| \leq a_1+a_2|x| \quad \forall x \in \R^n.
\end{align*}
In this case, testing the equation satisfied by $y_i$ with $y_i$, and a Gronwall argument yields
\begin{align*}
\left|y^k_i(t)\right|^2 \leq C_N \left(\left|y^i_0\right|^2+ \left\|y^k_i\right\|^2_{L^2} +
\left\|\F^\eps_{\theta_k}\left(y^k_i\right)\right\|^2_{L^2} \right)
\end{align*}
for some~$N$-dependent~$C_N>0$ and all~$t\in I$. Thus, the uniform boundedness of~$y^k_i$ in~$L^\infty(I;\R^n)$ follows again from the~$L^2$-estimates on~$y^k_i$ and~$\F^\eps_{\theta_k}(y^k_i)$ in~\eqref{eq:dimdepend}.
\end{itemize}
\end{rema}

The convergence result as $\eps\to 0^+$ of Theorem~\ref{thm:convoftraj} can be transferred to the finite training set setting as well.



\begin{prop}\label{prop:convergencefinite}
Let the regularisation parameters satisfy $\gamma_\eps
\|\theta_\eps\|^2_{\clR_\eps} = O(\eps)$. Further let~$\eps_k>0$ be a positive null sequence such that for each~$k \in \N$ there exists a solution~$(\bfy^k,\bfp^k,\theta_k) \in W_T^{2N} \times
\clR_\eps$ to~$(\clP^N_{\eps_k})$. If there is~${M_\infty}>0$ with~$\max_{i} \|y^k_i\|_{L^\infty} \leq M_\infty$ for all~$k \in \N$, then~$(\bfy^k,\bfp^k,
\mathbf{\F}^{\eps_k}_{\theta_k}(\bfy^k))$ admits at least one strong accumulation point~$(\mathbf{\bar{y}}, \mathbf{\bar{p}},
\mathbf{\bar{u}})$ in~$W^{2N}_T \times L^2(I;\R^m)^N$. Each such point satisfies
\begin{align*}
(\bar{y}_i, \bar{u}_i)\in \argmin \left(P^{y^i_0}_\beta\right), \quad i=1,\,\dots,\,N,
\end{align*}
as well as
\begin{align*}
\dot{\bar{y}}_i &=\bff(\bar{y}_i)+\bfg(\bar{y}_i)\bar{u}_i,~\bar{y}_i(0)=y^i_0
\\
-\dot{\bar{p}}_i&= D\bff(\bar{y}_i)^\top \bar{p}_i+ \left\lbrack D\bfg(\bar{y}_i)^\top\bar{u}_i\right\rbrack \bar{p}_i+\bfQ_1^\top
\bfQ_1(\bar{y}_i-y_d),~\bar{p}_i(T)= Q^\top_2 Q_2\left(\bar{y}_i(T)-y_d^T\right).
\end{align*}
\end{prop}

\begin{proof}
For every~$\eps_k$, with $k$ sufficiently large, denote by~$\theta_{\eps_k}\in \clR_{\eps_k}$ the corresponding parameters from Assumption~\ref{ass:approxsmoothness}, by~$\bfy_{\eps_{k}}$ the associated ensemble solution, see Theorem~\ref{thm:existenceneuralnetwork2}, and by~$\bfp_{\eps_{k}}$ the adjoint states. For abbreviation we set~$y^{\eps_k}_i \coloneqq \bfy_{\eps_k}(y^i_0)$ and $p^{\eps_k}_i \coloneqq \bfp_{\eps_k}(y^i_0)$. Then, by optimality, we have
\begin{align}\label{eq:aux11}
\sum^N_{i=1} \omega_i J\left(y^k_i, \F^{\eps_k}_{\theta_k}\left(y^k_i\right)\right) \leq
\sum^N_{i=1} \omega_i J_\eps\left(y^{\eps_k}_i, p^{\eps_k}_i,
\theta_{\eps_k}\right) + \frac{\gamma_{\eps_k}}{2} \left\|\theta_{\eps_k}\right\|^2_{\clR_\eps}.
\end{align}
As in the proof of Theorem~\ref{thm:convofobj} we see that the righthandside of this inequality converges to~$\sum^N_{i=1} \omega_i V^*(0,y^i_0)$ as~$k \rightarrow +\infty$. Thus it is bounded independently of~$k \in \N$. Similarly to Proposition~\ref{prop:existencefinite} we then conclude the existence of~$C_N >0$ depending on~$N$, but not on~$k$, such that
\begin{align*}
\max_{i} \left\|{\bfQ_1}y^k_i\right\|^2_{L^2}+ \max_{i}
\left\|\F^\eps_{\theta_k}\left(y^k_i\right)\right\|^2_{L^2} \leq C_N.
\end{align*}
Utilizing the state equation this can be improved to a~$k$-independent bound on the~$W_T$-norm of~$y^k_i$. By a Gronwall-type argument the same can be shown for the adjoint states~$p^k_i$. Now fix an arbitrary index~$i \in \{1,\,\dots,\,N\}$. Summarizing the previous observations we get the uniform boundedness of~$(y^k_i, p^k_i, \F^\eps_{\theta_k}(y^k_i))$ in~$W_T^2 \times L^2(I; \R^m) $ w.r.t. $k$, for each $i=1,\,\dots,\,N$. Each of its weak accumulation points~$(\bar{y}_i, \bar{p}_i, \bar{u}_i) \in W_T^2
\times L^2(I; \R^m)$ satisfies
\begin{align*}
\dot{\bar{y}}&=\bff(\bar{y})+\bfg(\bar{y})\bar{u},~
\bar{y}(0)=y_0.
\end{align*}
From this we conclude that
\begin{align*}
0\le \sum^N_{i=1} \omega_i V^*\left(0,y^i_0\right) \leq \sum^N_{i=1} \omega_i J(\bar y_i,\bar u_i) \le \lim_{k\,\to\,\infty} \sum^N_{i=1} \omega_i J\left(y^k_i, \F^{\eps_k}_{\theta_k}\left(y^k_i\right)\right)\le \sum^N_{i=1} \omega_i V^*\left(0,y^i_0\right),
\end{align*}
Since the second and third of the above inequalities also hold for each summand we conclude that~$\lim_{k\,\to\,\infty}J(y^k_i,
\F^{\eps_k}_{\theta_k}(y^k_i))\rightarrow J(\bar{y}_i, \bar{u}_i)$ as well as~$J(\bar{y}_i, \bar{u}_i)=V^*(0,y^i_0)$. Hence
\begin{align*}
(\bar{y}_i, \bar{u}_i)\in \argmin \left(P^{y^i_0}_\beta\right).
\end{align*}
The proof can now be concluded with minor adaptations to the proof of Theorem~\ref{thm:convoftraj}.
\end{proof}

A result analogous to that of Proposition~\ref{prop18} can also be obtained for Problem~\eqref{eq:learningprobfinite}. For the sake of brevity we do not present the details.


\subsection{The reduced objective functional}
In order to compute a solution to~\eqref{eq:learningprobfinite} we will rely on gradient-based optimization methods. For this purpose we introduce a~\emph{reduced objective functional} by eliminating the state and adjoint equations in~\eqref{eq:learningprobfinite}. Subsequently, we characterize the derivative of the reduced functional by means of adjoint techniques. To simplify the presentation we fix an arbitrary index~$i\in\{1,\,\dots,\,N\}$ in the following. Moreover, for abbreviation, we define the mapping
\begin{align*}
\bfA  \colon W_T \times \clR_\eps \to \clB \left(W_T;L^2\left(I;\R^n\right)\right),~A(y,\theta)=D\bff(y)^\top + \left\lbrack D\bfg(y)^\top\F^\eps_\theta(y)\right\rbrack.
\end{align*}
Using this notation, the adjoint equation in~\eqref{eq:learningprobfinite} can be expressed compactly as
\begin{align*}
-\dot{p}_i= \bfA (y_i,\theta)p_i+\bfQ_1^\top \bfQ_1(y_i-y_d),~p_i(T)=Q^\top_2 Q_2\left(y_i(T)-y_d^T\right).
\end{align*}
First, we argue the existence of~\emph{parameter-to-state operators} for the adjoint and the state equation.

\begin{lemm}\label{lem:existenceparametertostate}
Define~$G_i \colon W_T \times W_T \times \clR_\eps \to L^2(I; \R^n) \times L^2(I; \R^n) \times \R^n \times \R^n$ by
\begin{align*}
G_i(y,p, \theta)= \left(
\begin{array}{c}
\dot{y} -\bff(y)-\bfg(y)\F^\eps_\theta(y) \\
-\dot{p}- A(y,\theta)p-\bfQ_1^\top \bfQ_1(y-y_d)\\
y(0)-y^i_0 \\
p(T)- Q^\top_2 Q_2\left(y(T)-y_d^T\right)
\end{array}
\right).
\end{align*}
Let~$(\tilde{y},\tilde{p},\tilde{\theta})\in W_T \times W_T \times \clR_\eps$ satisfy~$G(\tilde{y},\tilde{p},\tilde{\theta})=0$. Then there exists a neighbourhood~$\clN_i(\tilde{y})\times \clN_i(\tilde{p}) \times \clN_i(\tilde{\theta})$ as well as~$\clC ^1$-mappings~$Y_i \colon \clN_i(\tilde{\theta})\to \clN_i(\tilde{y})\subset W_T$,~$P_i \colon \clN_i(\tilde{\theta})\to \clN_i(\tilde{p}) \subset W_T$ such that
\begin{align*}
G_i\left(Y_i(\theta),P_i(\theta),\theta\right)=0 \quad \forall \theta \in \clN(\tilde{\theta}).
\end{align*}
Given $y_i:=Y_i(\theta)$ and $~p_i:=P_i(\theta)$, the Fr\'{e}chet derivatives of~$Y_i$ and~$P_i$ at~$\theta \in \clN_i(\tilde{\theta})$, in direction~$\delta \theta \in \clR_\eps$, denoted by~$\delta Y_i \coloneqq Y'_i(\theta)(\delta \theta) $,~$\delta P_i \coloneqq P'_i(\theta)(\delta \theta)$ satisfy
\begin{multline*}
\dot{\delta Y_i}- \bfA (y_i,\theta)^\top \delta Y_i- \bfg(y_i) D_y \F^\eps_\theta(y_i)\delta Y\\
\begin{aligned}
&= \bfg(y_i) D_\theta \F^\eps_\theta(y_i)\delta \theta,
-\dot{\delta P_i}-\bfA (y_i,\theta)\delta P_i\\ &=\left\lbrack D_y\bfA (y_i,\theta)\delta Y_i \right\rbrack p_i+ \bfQ_1 \bfQ_1 \delta Y_i+\left\lbrack \partial_\theta \bfA (y_i,\theta) \delta \theta \right\rbrack p_i,\delta Y_i(0)\\
&=0,~\delta P_i(T)= Q_2^\top Q_2 \delta Y_i(T).
\end{aligned}
\end{multline*}
\end{lemm}
\begin{proof}
This is a direct consequence of the implicit function theorem applied to~$G$ noting that the directional derivatives satisfy
\begin{align*}
\begin{pmatrix}
\partial_y G_i(y,p,\theta) & \partial_p G_i(y,p,\theta)
\end{pmatrix}
\begin{pmatrix}
\delta Y \\ \delta P
\end{pmatrix} =- \partial_\theta G_i(y,p,\theta) \delta \theta.
\end{align*}
\end{proof}

Now consider an admissible point~$(\tilde{\bfy}, \tilde{\bfp}, \theta) \in W^{2N}_T \times \clR_\eps$ for~\eqref{eq:learningprobfinite}. For every~$i=1,\,\dots,\,N$, let~$\clN_i(\tilde{\theta})$ and~$Y_i,P_i$ denote the corresponding neighbourhoods and operators from Lemma~\ref{lem:existenceparametertostate}. Setting~$\clN(\tilde{\theta})= \bigcap^N_{i=1} \clN_i(\tilde{\theta})$ define the reduced objective functional
\begin{align}\label{eq:reducedobj}
\clJ_N \colon \clN(\tilde{\theta}) \to [0,+\infty),~\clJ_N(\theta)= \sum^N_{i=i} \omega_i J_\eps\left(Y_i(\theta),P_i(\theta), \theta\right)+ \frac{\gamma_\eps}{2}\|\theta\|^2_{\clR_\eps},
\end{align}
and set
\begin{align*}
\Phi_i(t)= \int^t_0 \left(V^\eps_\theta(s,y_i(s))-J_s(y_i,u)\right)~\de s.
\end{align*}


\begin{prop}\label{prop:direcderiv}
The functional~$\clJ_N$ from~\eqref{eq:reducedobj} is at least of class~$\clC ^1$ on~$\clN(\tilde{\theta})$. Given~$ \theta \in \clN(\tilde{\theta})$, set~$y_i \coloneqq Y_i(\theta)$,~$p_i \coloneqq P_i(\theta)$ as well as~$\delta Y_i \coloneqq Y'_i(\theta)(\delta \theta)$,~$ \delta P_i \coloneqq P'_i(\theta)(\delta \theta)$. The directional derivative of~$\clJ_N$ at~$\theta$ in the direction of~$\delta \theta \in \clR_\eps$ is given by
\begin{align*}
\clJ'_N(\theta)(\delta \theta)= \sum^N_{i=1} \omega_i \left((\widehat{y}_i, \delta Y_i)_{L^2}+\left(\widehat{y}^T_i, \delta Y_i(T)\right)_{\R^n}+\left(\widehat{p}_i, \delta P_i\right)_{L^2}+\left(\widehat{\theta}_i,\delta \theta\right)_{\clR_\eps}\right)+\gamma_\eps (\theta,\delta \theta)_{\clR_\eps}
\end{align*}
with
\begin{multline*}
\widehat{y}_i=(1-\gamma_1\Phi_i) \bfQ_1 \bfQ_1 (y_i-y_d)+ \beta(1-\gamma_1\Phi_i) D_y \F^\eps_\theta(y_i)^\top\F^\eps_\theta(y_i)\\+\gamma_1 \left(\clV^\eps_\theta(t,y_i)-J_\bullet\left(y_i,\F^\eps_\theta(y_i)\right)\right)\partial_y \clV^\eps_\theta(y_i)+\gamma_2 D_{yy}\clV^\eps_\theta(y_i)\left(\partial_{y}\clV^\eps_\theta(y_i)-p_i\right),
\end{multline*}
and
\begin{align*}
\widehat{y}^T_i=\alpha\left(1-\gamma_1\Phi_i(0)\right) Q_2 Q_2 \left(y_i(T)-y^T_d\right),
\end{align*}
as well as
\begin{align*}
\widehat{p}_i=\gamma_2\left(p_i-\partial_y \clV^\eps_\theta(y_i)\right),
\end{align*}
and
\begin{multline*}
\widehat{\theta}_i = \gamma_1 \int^T_0 D_\theta V^\eps_\theta(t,y_i(t))^\top\left(V^\eps_\theta(t,y_i(t))-J_t\left(y_i,\F^\eps_\theta(y_i)\right)\right)~\de  t \\
+\int^T_0 \!\left\lbrack \beta \left(1\!-\!\gamma_1\Phi_i(t)\right) D_\theta F^\eps_\theta(t,y_i(t))^\top F^\eps_\theta(t,y_i(t))\!+\!\gamma_2 D_{y \theta} V^\eps_\theta(t,y_i(t))^\top\left(\partial_{y}{V}^\eps_\theta(t,y_i(t)\right)\!-\!p_i(t)) \right\rbrack~\de  t.
\end{multline*}
\end{prop}

\begin{proof}
The regularity of~$\clJ_N$ follows immediately from Lemma~\ref{lem:existenceparametertostate} and the chain rule. In order to compute the directional derivative we abbreviate
\begin{align*}
F_1(y,u,\theta)=\frac{\gamma_1}{2}\int^T_0 \left|V^\eps_\theta(t,y(t))-J_t(y,u)\right|^2\,\de  t,\\ F_2(y,p,\theta)=\frac{\gamma_2}{2}\int^T_0 \left|\partial_y V^\eps_\theta(t,y(t))-p(t)\right|^2~\de  t
\end{align*}
in the following. Thus we have
\begin{align*}
J_\eps(Y_i(\theta),P_i(\theta),\theta)&= J\left(Y_i(\theta),\F^\eps_\theta(Y_i(\theta))\right)+F_1\left(Y_i(\theta),\F^\eps_\theta(Y_i(\theta)),\theta\right)+F_2\left(Y_i(\theta),P_i(\theta),\theta\right)
\\
&=G_1(\theta)+G_2(\theta)+G_3(\theta).
\end{align*}
We readily verify
\begin{multline*}
G'_1(\theta)(\delta \theta)=\left(\bfQ_1 \bfQ_1 (y_i-y_d),\delta Y_i\right)_{L^2}+ \beta \left(D_y \F^\eps_\theta(y_i)^\top\F^\eps_\theta(y_i),\delta Y_i\right)_{L^2}\\
+ \beta \left(D_\theta \F^\eps_\theta(y_i)^\top\F^\eps_\theta(y_i),\delta \theta\right)_{\clR_\eps}+ \left(Q_2 Q_2 \left(y_i(T)-y^T_d\right),\delta Y_i(T)\right)_{\R^n}.
\end{multline*}
Recalling the definition of $\Phi_i$ we get
\begin{align*}
G'_2(\theta)& (\delta \theta)=\gamma_1 \left(E_1+E_2+E_3+E_4\right),
\end{align*}
where
\begin{align*}
E_1 &= \left(\left(\clV^\eps_\theta(y_i)\!-\!J_\bullet\left(y_i,\F^\eps_\theta(y_i)\right)\right) \partial_y \clV^\eps_\theta(y_i), \delta Y_i\right)_{L^2} +\left(D_\theta \clV^\eps_\theta(y_i)^\top\left(\clV^\eps_\theta(y_i)\!-\!J_\bullet\left(y_i,\F^\eps_\theta(y_i)\right)\right), \delta \theta\right)_{\clR_\eps} \\
&= \left(\left(\clV^\eps_\theta(y_i)\!-\!J_\bullet\left(y_i,\F^\eps_\theta(y_i)\right)\right) \partial_y \clV^\eps_\theta(t,y_i), \delta Y_i\right)_{L^2} \\
&\qquad +\left(\int^T_0 D_\theta V^\eps_\theta\left(t,y_i(t)\right)^\top\left(V^\eps_\theta(t,y_i(t))-J_t\left(y_i,\F^\eps_\theta(y_i)\right)\right)~\de  t, \delta \theta \right)_{\clR_\eps},
\end{align*}
\begin{align*}
E_2 &=- \int^T_0 \left(V^\eps_\theta(t,y(t))-J_t(y,u)\right)\\
&\qquad\left(\int^T_t \left(Q_1 Q_1(y(s)-y_d(s)),\delta y(s)\right)~\de  s+\left(Q_2 Q_2\left(y(T)-y^T_d\right),\delta y(T)\right)_{\R^n} \right)~\de  t \\
&=-\left(\Phi_i Q_1 Q_1(y-y_d),\delta y\right)_{L^2}- \Phi_i(0) \left(Q_2 Q_2(y(T)-y^T_d\right),\delta y(T))_{\R^n},
\end{align*}
as well as
\begin{align*}
E_3&=- \int^T_0 \left(V^\eps_\theta(t,y(t))-J_t(y,u)\right) \left(\beta \int^T_t\left(D_y F^\eps_\theta(s,y_i(s))^\top F^\eps_\theta(s,y_i(s)),\delta Y_i(s)\right)_{\R^n}~\de  s\right)~\de  t \\ &=- \beta\left(\Phi_i D_y \F^\eps_\theta(y_i)^\top \F^\eps_\theta(y_i), \delta Y_i\right)_{L^2},
\end{align*}
and
\begin{align*}
E_4&=- \int^T_0 \left(V^\eps_\theta(t,y(t))-J_t(y,u)\right) \left(\beta \int^T_t \left(D_\theta F^\eps_\theta\left(s,y_i(s)\right)^\top F^\eps_\theta(s,y_i(s)),\delta \theta\right)_{\clR_\eps}~\de  s\right)~\de  t \\ &=- \beta \left(\int^T_0 \Phi_i(t) D_\theta F^\eps_\theta(t,y_i(t))^\top F^\eps_\theta(t,y_i(t))~\de  t, \delta \theta \right)_{\clR_\eps},
\end{align*}
by means of partial integration. Finally we calculate
\begin{multline*}
G'_3(\theta)(\delta \theta)= \gamma_2 \left(D_{yy}\clV^\eps_\theta(y_i)\left(\partial_{y}\clV^\eps_\theta(y_i)-p_i\right), \delta Y_i\right)_{L^2}- \gamma_2\left(\partial_{y}\clV^\eps_\theta(y_i)-p_i, \delta P_i\right)_{L^2} \\ +\gamma_2\left(\int^T_0 D_{y \theta} V^\eps_\theta(t,y_i(t))^\top\left(\partial_{y}{V}^\eps_\theta(t,y_i(t))-p_i(t)\right)~\de  t, \delta \theta \right)_{\clR_\eps}.
\end{multline*}
Summarizing the previous observations, we arrive at the claimed characterization.
\end{proof}

Applying a gradient method to~\eqref{eq:learningprobfinite} requires the computation of the gradient~$\nabla \clJ_{N}(\theta)\in \clR_\eps$ which satisfies
\begin{align*}
\clJ'_N(\theta)(\delta \theta)=\left(\nabla \clJ_N(\theta), \delta \theta\right)_{\clR_\eps} \quad \forall \delta \theta \in \clR_\eps.
\end{align*}
This can be done by computing~$\clJ'_N(\theta)(e_j)$ for the canonical basis~$\{e_j\}^{N_\eps}_{j=1}\subset \clR_\eps$. However, such reasoning leads to the necessity to solve~$2\dim(\clR_\eps) N$ additional ODEs in order to compute the sensitivities~$Y'_i(\theta)(e_j)$ and~$P'_i(\theta)(e_j)$, respectively. Introducing suitable costate equations, this can be reduced to~$2N$ additional equation solves.

\begin{lemm}\label{lem:adjoint}
Let~$\widehat{y}_i,\widehat{y}^T_i, \widehat{p}_i$ as well as~$\delta Y_i, \delta P_i$ be defined as in Proposition~\ref{prop:direcderiv}. Then there holds
\begin{align*}
\left(\widehat{y}_i,\delta Y_i\right)_{L^2}+ \left(\widehat{y}^T_i,\delta Y_i (T)\right)_{\R^n}+\left(\widehat{p}_i,\delta P\right)_{L^2}=\left(D_\theta \F^\eps_\theta(y_i)^\top\left(\bfg(y_i)^\top \zeta_i+\left\lbrack D\bfg(y_i) \kappa_i \right\rbrack^\top p_i\right), \delta \theta\right)_{\clR_\eps}
\end{align*}
where~$\zeta_i, \kappa_i \in W_T$ satisfy
\begin{align*}
-\dot{\zeta}_i&=\bfA (y_i,\theta)\zeta+ D_y \F^\eps_\theta(y_i)^\top \bfg(y_i)^\top\zeta_i + \left\lbrack D_y \bfA (y_i,\theta)^\top p\right\rbrack\kappa_i+\bfQ_1^\top \bfQ_1\kappa_i+\widehat{y}_i \\
\dot{\kappa}_i&=\bfA (y_i,\theta)^\top \kappa_i+\widehat{p}_i,\\
\zeta_i(T)&= Q_2^\top Q_2 \kappa (T)+\widehat{y}^T_i,~\kappa_i(0)=0.
\end{align*}
\end{lemm}

\begin{proof}
For the sake of readability, we drop the subscript~$i$ in the following. By partial integration and Lemma~\ref{lem:existenceparametertostate} we obtain
\begin{align*}
\left(\widehat{p},\delta P\right)_{L^2}&=\left(\dot{\kappa}-\bfA (y,\theta)^\top \kappa,\delta P\right)=\left(-\dot{\delta P}-\bfA (y,\theta) \delta P,\kappa\right)+\left(Q_2^\top Q_2 \delta Y(T),\kappa(T)\right)_{\R^n}\\
&=\left(\left\lbrack D_y\bfA (y,\theta)\delta Y \right\rbrack p+ \bfQ_1^\top \bfQ_1 \delta Y+\left\lbrack \partial_\theta \bfA (y,\theta) \delta \theta \right\rbrack p,\kappa\right) +\left(\delta Y(T),\zeta(T)- \hat y_i^T\right)_{\R^n}
\end{align*}
and
\begin{multline*}
(\hat y_i,\delta Y)_{L^2}+\left(\hat y_i^T,\delta Y (T)\right)_{\R^n}\\
\begin{aligned}
&=\left(-\dot{\zeta}\!-\!\bfA (y,\theta)\zeta\!-\! D_y \F^\eps_\theta(y)^\top \bfg(y)^\top\zeta- \left\lbrack D_y \bfA (y,\theta)^\top p\right\rbrack\kappa-\bfQ_1^\top \bfQ_1\kappa, \delta Y\right)_{L^2}+ \left(\hat y^T_i,\delta Y (T)\right)_{\R^n}\\
&=\left(\dot{\delta Y}-\bfA (y,\theta)^\top\delta Y- \bfg(y) D_y \F^\eps_\theta(y) \delta Y,\zeta\right)_{L^2}\\
&\quad -\left(\left\lbrack D_y \bfA (y,\theta)^\top p\right\rbrack\kappa+\bfQ_1^\top \bfQ_1\kappa, \delta Y\right)_{L^2} - \left(\delta Y (T), \zeta(T)\right)_{\R^n} \\
&=\left(\bfg(y) D_\theta \F^\eps_\theta(y)\delta \theta,\zeta\right)_{L^2} -\left(\left\lbrack D_y \bfA (y,\theta)^\top p\right\rbrack\kappa+\bfQ_1^\top \bfQ_1\kappa, \delta Y\right)_{L^2} - \left(\delta Y (T), \zeta(T)- \hat y_i^T\right)_{\R^n}.
\end{aligned}
\end{multline*}
Adding both equations finally yields
\begin{align*}
\left(\hat y_i,\delta Y\right)_{L^2}+ \left(\hat y_i^T,\delta Y (T)\right)_{\R^n}+\left(p_1,\delta P\right)_{L^2}&=\left(\bfg(y) D_\theta \F^\eps_\theta(y)\delta \theta,\zeta\right)_{L^2} +\left(\left\lbrack \partial_\theta \bfA (y,\theta) \delta \theta \right\rbrack p,\kappa\right)_{L^2}
\\
&=\left(D_\theta\F^\eps_\theta(y)^\top\left(\bfg(y)^\top \zeta+\lbrack D\bfg(y) \kappa \rbrack^\top p\right), \delta \theta\right)_{\clR_\eps}
\end{align*}
which ends the proof of Lemma~\ref{lem:adjoint}.
\end{proof}

We arrive at the following characterization of the gradient~$\nabla \clJ_N(\theta)$.

\begin{theo}\label{thm:gradient}
Let~$y_i,p_i,\zeta_i, \kappa_i \in W_T,~\widehat{\theta}_i \in \clR_\eps$ be defined as in~Proposition~\ref{prop:direcderiv} and Lemma~\ref{lem:adjoint}. The gradient of~$\clJ_N$ at~$\theta$ is given by
\begin{align*}
\nabla \clJ_N(\theta)= \sum^N_{i=1} \omega_i \left(D_\theta \F^\eps_\theta(y_i)^\top\left(\bfg(y_i)^\top \zeta_i+\left\lbrack D\bfg(y_i) \kappa_i \right\rbrack^\top p_i\right) +\widehat{\theta}_i\right)+ \gamma_\eps \theta.
\end{align*}
\end{theo}
%~$\mu$ is a regular Borel measure and $\bfy_{\theta^*_{\eps,T}}\in \clC _b(\supp \mu;W_T)$ then all of the results in this section can be proven replacing

\section{Numerical example}\label{sec:numericalexampl}
We finish this paper by applying the proposed learning approach to one particular instance of Problem~\eqref{def:openloopproblem}. Setting~$I=(0,T)$ and~$\Omega=(0,2\pi)$, we consider the parabolic bilinear optimal control problem
\begin{align*}
\min_{\clY\,\in\,L^2(I \times \Omega),\,u\,\in\,L^2\left(I;\R^3\right)} \left \lbrack \frac{1}{2} \int_I \left\|\clY (t)-\clY _d(t)\right\|^2_{L^2(\Omega)}+ \frac{\beta}{2}|u(t)|^2_{\R^3}~\de t \right\rbrack+\frac{\alpha}{2} \left\|\clY (T)-\clY _d(T)\right\|^2_{L^2(\Omega)}
\end{align*}
subject to
\begin{align}\label{eq:bilinearPDE}
\partial_t \clY -\bigtriangleup \clY + \left(u_1 \chi_1 + u_2 \chi_2+ u_3 \chi_3 \right) \clY = 0,
\end{align}
as well as
\begin{align*}
\clY (t,x)=0 \quad \text{on}~I \times \partial \Omega,~\clY (0,x)= \clY _0(x) \quad \text{on}~\Omega.
\end{align*}
Here~$\alpha>0, \beta>0$, and~$\clY _d$ denotes a given desired state. The dynamics of this infinite-dimensional system can be influenced by choosing a time-dependent three-dimensional control input~$u \in L^2(I; \R^3) $ which acts on the subdomains~$\Omega_1=(0.5,1)$,~$\Omega_2=(2,2.5)$ and $\Omega_3=(4,4.5)$, respectively. The associated characteristic functions are denoted by~$\chi_i$,~$i=1,\,\dots,\,3$.

In order to fit this problem into the setting of the current manuscript, let~$\{\lambda_i, \varphi_i\} \in \R_+ \times L^2(\Omega)$ denote the first $n \in \N$ normalized eigenpairs of the Dirichlet Laplacian on~$\Omega$. Approximating the state dynamics~$\clY $ as well as the desired state by
\begin{align*}
\clY (t,x) \approx \sum^n_{i=1} Y_i(t) \varphi_i (t),~\clY _d(t,x) \approx \sum^n_{i=1} Y^i_d(t) \varphi_i (t),~
\end{align*}
we end up with
\begin{align}\label{def:discreteexample}
\min_{{Y}\,\in\,L^2\left(I;\R^{10}\right),\,u\,\in \,L^2\left(I;\R^3\right)} \left \lbrack \frac{1}{2} \int_I \left|Y(t)-Y_d(t)\right|^2+ \frac{\beta}{2}|u(t)|^2_{\R^3}~\de t+ \frac{\alpha}{2}\left|Y(T)-Y_d(T)\right|^2_{\R^{10}} \right \rbrack
\end{align}
subject to
\begin{align*}
\dot{Y}(t)+ AY(t)+ \sum^3_{i=1} u_i M_i Y(t) =0,~Y(0)=Y_0.
\end{align*}
where~$(Y_0)_i= (\clY _0, \varphi_i)_{L^2}$,~$i=1,\,\dots,\,n$, and the symmetric matrices~$A,M_i \in \R^{n\times n}$ are given by
\begin{align*}
A_{jk}=
\begin{cases} 0 & j \neq k \\
\lambda_j & \text{else}
\end{cases},~(M_i)_{jk}=\int_\Omega \phi_j \phi_k \chi_i(x)~\de x, \quad i=1,2,3,~j,k= 1,\dots, n.
\end{align*}
Note that by choosing a spectral technique to approximate the infinite dimensional system~\eqref{eq:bilinearPDE} by a finite dimensional one, we have done justice to the fact that grid based techniques would rapidly lead to systems of dimensions which are challenging or even impossible. Special techniques, as for instance tensor train methods, would then become essential.


\subsection{Learning \& validation setup}
In the following, we determine an approximate optimal feedback law for~\eqref{def:discreteexample} by applying the learning approach detailed in Section~\ref{sec:learnfeedback}. The parametrized model~$V^\eps_\theta$ for the value function is given by realizations of residual networks, as described in Section~\ref{subsec:residual}, with~$L_\eps=2$ layers,~$\arch(\theta)=(11,60,1)$ and activation function~$\sigma$ given by
\begin{align*}
\sigma(x)= \sin(x)+\cos(x).
\end{align*}
This yields a total of~$1440$ trainable parameters. We emphasize that the architecture as well as the activation function were chosen based on numerical testing. In particular, the present tests should not be mistaken as a \emph{quantitative} survey but as a \emph{proof of concept} which highlights the potential of learned feedbacks for optimal control and puts a focus on the role played by the penalty parameters~$\gamma_1$ and~$\gamma_2$.

Given a fixed reference vector~$\bar{Y}_0 $, we randomly generate a set~$\bfy_0$ of~$130$ initial conditions by sampling uniformly from the closure of~$B_1(\bar{Y}_0)$, Subsequently, these are split into a training set~$\bfy^t_0$ of $N=30$ initial conditions, which is used in the learning problem~\eqref{eq:learningprobfinite} together with uniform weights~$w_j=1/N$, and a validation set~$\bfy^v_0= \bfy_0 \setminus \bfy^T_0$ which we later utilize to assess the performance of the obtained feedback.

In order to obtain a candidate for the optimal network parameters~$\theta^*_\eps$, a Barzilai-Borwein method~\cite{BaBo19}, is applied to the learning problems~\eqref{eq:learningprobfinite}, based on the reduced objective functional introduced in~\eqref{eq:reducedobj} as well as the characterization of its gradient in Theorem~\ref{thm:gradient}. For every~$Y_0 \in \bfy^t_0$, this approach entails the computation of the state~$Y \coloneqq Y_\theta (Y_0)$ and the adjoint state~$P \coloneqq P_\theta (Y_0)$ which satisfy
\begin{equation}\label{def:exampleneuralequations}
\begin{split}
\dot{Y} (t)+ \left(A+ \sum^3_{i=1} {F}^\eps_\theta(t,Y(t))_i M_i \right) Y(t) &=0,~Y(0)=Y_0 \\
-\dot{P} (t) + \left(A+ \sum^3_{i=1} {F}^\eps_\theta(t,Y(t))_i M_i \right) P(t) &= Y(t)-Y_D(t),~P(T)= Y(T)-Y_D(T)
\end{split}
\end{equation}
as well as the costates~$K \coloneqq K_\theta(Y_0)$ and~$Z \coloneqq Z_\theta(Y_0)$ with
\begin{align*}
\dot{K}(t)+ \left(A+ \sum^3_{i=1} {F}^\eps_\theta(t,Y(t))_i M_i \right) K(t)= \widehat{P}(t)
\end{align*}
and
\begin{multline*}
- \dot{Z} (t)+ \left(A+ \sum^3_{i=1} {F}^\eps_\theta(t,Y(t))_i M_i + D_y {F}^\eps_\theta(t,Y(t))^\top
\begin{pmatrix}
Y_j(t)^\top M_1 \\
Y_j(t)^\top M_2 \\
Y_j(t)^\top M_3
\end{pmatrix} \right) Z(t) \\
=- D_y {F}^\eps_\theta(t,Y(t))^\top
\begin{pmatrix}
Y(t)^\top M_1 \\
Y(t)^\top M_2 \\
Y(t)^\top M_3
\end{pmatrix} Z (t) +K (t)+ \widehat{Y} (t)
\end{multline*}
equipped with the boundary conditions
\begin{align*}
K(0)=0, \quad Z(T)= \alpha K(T)+ \widehat{Y}^T_j
\end{align*}
where~$\widehat{Y},\widehat{Y}^T$ and~$\widehat{P}$ are defined in analogy to Proposition~\ref{prop:direcderiv}. Note that this system is not fully coupled, i.e. in practice, we first solve the nonlinear closed-loop equation using a Radau time-stepping scheme and then, successively treat the adjoint and costate equations by an implicit Euler method. This can be done in parallel for various initial conditions to achieve additional speed-up. Moreover, the adjoint state~$P$ and costate~$K$ only need to be computed if~$\gamma_2 >0$. The gradient of the reduced objective functional~$\clJ_N$ in~\eqref{eq:learningprobfinite} at an admissible~$\theta$ is then obtained as
\begin{align*}
\frac{1}{30} \sum_{Y \in \bfy^t_0} \int_I \left(D_\theta F^\eps_\theta(t,Y_\theta (Y_0)(t))^\top \left(B^\theta_Y (t) Z_\theta (Y_0) (t)+ B^\theta_K(t) P_\theta (Y_0) (t) \right)~\de t +\widehat{\theta}(Y_0)\right).
\end{align*}
where we set
\begin{align*}
B^\theta_Y (t) \coloneqq
\begin{pmatrix}
Y_\theta (Y_0)(t)^\top M_1 \\
Y_\theta (Y_0)(t)^\top M_2 \\
Y_\theta (Y_0)(t)^\top M_3
\end{pmatrix}, \quad B^\theta_K(t) \coloneqq
\begin{pmatrix}
K_\theta (Y_0)(t)^\top M_1 \\
K_\theta (Y_0)(t)^\top M_2 \\
K_\theta (Y_0)(t)^\top M_3
\end{pmatrix},
\end{align*}
integration has to be understood componentwise and~$\widehat{\theta}(Y_0)$ is as in Proposition~\ref{prop:direcderiv}.


Once the network is determined, we compute the state~$Y_{\theta} (Y_0)$ and adjoint~$P_\theta (Y_0)$ for every~$Y_0 \in \bfy_0$ from~\eqref{def:exampleneuralequations} and set~$U_\theta(Y_0)\coloneqq\F^\eps_{\bar{\theta}}(Y_\theta(Y_0))$. Subsequently we determine a stationary point~$(\bar{Y}(Y_0), \bar{U}(Y_0))$ of~\eqref{def:discreteexample},~$Y_0 \in \bfy_0$, by applying a Barzilai--Borwein gradient method to its control-reduced formulation. The associated adjoint state is denoted by $\bar{P}(Y_0)$. At this point, it should be stressed that both, the open loop as well as the feedback learning problem, are nonconvex. As a consequence, we cannot ensure global optimality of the computed stationary points and, in particular, both methods might provide different results. For the present example, open loop and learned feedback controls are comparable. Moreover, for every~$Y_0 \in \bfy_0$, we have~$J(\bar Y(Y_0), \bar U(Y_0)) \geq J(Y_\theta (Y_0), U_\theta (Y_0))$. In order to assess the performance of open loop and feedback controls, let~$Y^{ad}_{0}\subset \bfy_0$ be either~$Y^{ad}_{0}=\bfy^t_0$ or~$Y^{ad}_{0}=\bfy^v_0$ and consider the relative difference between the averaged objective functional values:
\begin{align*}
\Errr_{\clJ} \coloneqq \frac{\sum_{Y_0\,\in\,Y_{ad}} J\left(Y_\theta (Y_0), U_\theta (Y_0)\right)-\sum_{Y_0\,\in\,Y_{ad}} J\left(\bar Y (Y_0),\bar U (Y_0)\right) }{\sum_{Y_0\,\in\,Y_{ad}} J\left(\bar Y (Y_0),\bar U (Y_0)\right)}
\end{align*}
as well as the associated normalized mean squared error of~$J(Y_\theta (\cdot), U_\theta (\cdot))$:
\begin{align*}
\Errr_{J} \coloneqq \frac{\sum_{Y_0\,\in\,Y_{ad}} \left(J\left(Y_\theta (Y_0), U_\theta (Y_0)\right)- J\left(\bar Y (Y_0),\bar U (Y_0)\right)\right)^2}{\sum_{Y_0\,\in\,Y_{ad}} J\left(\bar Y (Y_0),\bar U (Y_0)\right)^2}.
\end{align*}
The normalized mean-squared errors of the state, $\Errr_{Y}$, adjoint, $\Errr_{P}$, and of the control,~$\Errr_{U}$, are defined analogously. Moreover, to quantify the influence of the penalty parameters~$\gamma_1$ and~$\gamma_2$, we define
\begin{align*}
\Errr_{V} \coloneqq \frac{\sum_{Y_0\,\in\,Y_{ad}} \int_I \left|V^\eps_\theta \left(t,Y_\theta(Y_0)(t)\right)- J_t\left(Y_\theta(Y_0),U_\theta(Y_0)(t)\right)\right|^2~\de t}{\sum_{Y_0\,\in\,Y_{ad}} \int_I \left| J_t\left(Y_\theta(Y_0),U_\theta(Y_0)(t)\right)\right|^2~\de t}.
\end{align*}
as well as
\begin{align*}
\Errr_{\partial V} \coloneqq \frac{\sum_{Y_0\,\in\,Y_{ad}} \int_I \left|\partial_y V^\eps_\theta \left(t,Y_\theta(Y_0)(t)\right)- P_\theta (Y_0)(t)\right|^2~\de t }{\sum_{Y_0\,\in\,Y_{ad}} \int_I\left| P_\theta (Y_0)(t)\right|^2~\de t}.
\end{align*}
For $Y^{ad}_0= Y^t_0$, these terms correspond to the relative sizes of the additional penalties in~\eqref{eq:learningprobfinite}. Finally, we also want to compare~$V^\eps_\theta$ with the optimal value function~$V^*$. Of course,~$V^*$ can neither be given analytically nor can it be computed exactly. As a remedy, we recall that if~$V^*$ is sufficiently regular and~$(\bar{Y}(Y_0), \bar{U}(Y_0))$ is a minimizing pair of~\eqref{def:discreteexample} with adjoint state~$\bar{P}(Y_0)$, we have
\begin{align*}
V^*\left(t, \bar{Y}(Y_0)(t)\right)= J_t\left(\bar{Y}(Y_0),\bar{U}(Y_0)\right) \quad \text{as well as} \quad \partial_y V^*\left(t, \bar{Y}(Y_0)(t)\right)= \bar{P}(Y_0)(t)
\end{align*}
for all~$t \in I$. As a consequence, setting
\begin{align*}
d \left(V^*, V^\eps_\theta\right)=\frac{\sum_{Y_0\, \in\,Y_{ad}} \int^T_0 \left| V^\eps_\theta\left(t, \bar{Y}(Y_0)(t)\right)- J_t\left(\bar{Y}(Y_0),\bar{U}(Y_0)\right) \right|^2~\de t}{\sum_{Y_0\,\in\, Y_{ad}} \int^T_0 \left|J_t\left(\bar{Y}(Y_0),\bar{U}(Y_0)\right)\right|^2~\de t}.
\end{align*}
as well as
\begin{align*}
d \left(\partial V^*, \partial V^\eps_\theta\right)=\frac{\sum_{Y_0\,\in\,Y_{ad}} \int^T_0 \left| \partial_y V^\eps_\theta\left(t, \bar{Y}(Y_0)(t)\right)- \bar{P}(Y_0)(t)\right|^2~\de t}{\sum_{Y_0\,\in\,Y_{ad}} \int^T_0\left|\bar{P}(Y_0)(t)\right|^2~\de t}.
\end{align*}
provides a suitable ``distance'' for the comparison of~$V^*$ and~$V^\theta_\eps$.

\subsection{Validation results} 

As a concrete example, we set~$T=2,~\beta=0.01,~\alpha=0.25 $ and~$\clY _d(t,x)= x^2/10$, i.e., we try to steer the system towards a parabola. Note that there is no control input~$u \in L^2(I;\R^3)$ such that the corresponding solution~$\clY $ of the PDE~\eqref{eq:bilinearPDE} satisfies~$\clY (t)= \clY _d$. The parabolic binlinear control problem is approximated using~$n=10$ eigenfunctions. All computations were carried out in Matlab 2019 on a notebook with~$32$ GB RAM and an Intel\textregistered Core\texttrademark~i7-10870H CPU@2.20 GHz.


In order to compute an approximately optimal feedback law for this problem, we solve~\eqref{eq:learningprobfinite} for various penalty parameter configurations~$\gamma_1, \gamma_2 \in \{0,0.1,1\}$. The resulting normalized errors can be found in Table~\ref{tab:val1}, for~$Y^{ad}_0=\bfy^t_0$, and Table~\ref{tab:val2}, for~$Y^{ad}_0=\bfy^v_0$. Comparing their individual entries, we observe that there is (almost) no difference in performance between the training and the validation sets. This means that, while the utilized networks are rather simple and only comprise a small number of trainable parameters, the corresponding learned feedback controls generalize well to initial conditions which are not contained in the training set.


Indeed, on the one hand \emph{all} computed networks provide feedback controls which perform similarly to their open loop counterparts. This is manifested in very small averaged errors for the objective functional, i.e.~$\Errr_{\clJ}$ and $\Errr_{{J}}$, the states and adjoint states,~$\Errr_{{Y}}$ and $\Errr_{{P}}$, as well as the controls,~$\Errr_{{U}}$. These start to (slowly) deteriorate as~$\gamma_1$ and/or~$\gamma_2$ grow. However, cf. the explanation in Section~\ref{sec:learnfeedback}, this is expected: for~$\gamma_1>0$ and/or~$\gamma_2 >0$, the learned feedback has to strike a balance between minimizing~$J(Y_\theta (\cdot), U_\theta(\cdot))$ and keeping the penalty terms small, hence the slightly larger error.

On the other hand, the picture looks different once we consider the errors associated to the approximation of the value function, i.e.,~$\Errr_{V}$,~$\Errr_{\partial V}$ as well as $d(\partial V^*, \partial V^\eps_\theta)$ and $d (\partial V^*, \partial V^\eps_\theta)$. Here $\gamma_1>0$ and/or~$\gamma_2 >0$ have a significant influence on $d(V^*, V^\eps_\theta)$ and~$d(\partial V^*, \partial V^\eps_\theta)$ while the other normalized mean squared errors remain relatively small. Moreover, we have~$\Errr_{V} \approx d (V^*, V^\eps_\theta) $ and~$\Errr_{\partial V} \approx d (\partial V^*, \partial V^\eps_\theta) $ on the test as well as on the validation set. Hence, large values for these terms are a reliable indicator for structural differences between~$V^\eps_\theta$ and~$V^*$ and/or~$\partial_y V^\eps_\theta$ and~$ \partial_y V^*$, respectively.

Now, while~$\gamma_1=\gamma_2=0$ provides a very good approximation to the open loop optimal control, it performs the worst in terms of approximating the optimal value function and its derivative. This is related to two observations. First, in this case, the learning problem~\eqref{eq:learningprobfinite} only depends on the derivative~$\partial_y V^\theta_\eps$ but \emph{not} on the value function~$V^\eps_\theta$. Since primitives are not unique, approximating~$V^*$ by~$V^\eps_\theta$ is unlikely. Second, due to the absence of~$V^\eps_\theta$ in the problem, some of the parameters in the model are not trainable. In fact, for $\gamma_1=\gamma_2=0$, there holds~$\partial_{W_{12}} \clJ_N(\theta)=0$ for every admissible~$\theta$.

Once we increase~$\gamma_1$ and~$\gamma_2$, this is no longer the case. Hence, we observe rapid decrease for $d(V^*, V^\eps_\theta)$ and $d(\partial V^*,\partial V^\eps_\theta)$. Most remarkably, the improvement for both is, to some extend, already visible for~$\gamma_1>0$ and~$\gamma_2=0$. In this setting, applying the gradient method neither requires computing the adjoint state~$P$ nor the costate~$K$ which limits the cost of every gradient step to~$2N=60$ ODE solve. Quite the contrary, increasing~$\gamma_2 >0$ but keeping~$\gamma_1=0$ fixed, there is \emph{no} improvement for $d(\ V^*, V^\eps_\theta)$. This further backs up our reasoning given for the case of~$\gamma_1=\gamma_2=0$.

Consequently, the computed results indicate that the best balance between finding an optimal control and approximating the value function is achieved by a careful choice of~$\gamma_1,\gamma_2 >0$. Moreover, they highlight two important points: First, the presented learning approach indeed allows to compute semiglobal optimal feedback laws~$F^\eps_\theta$ for higher dimensional problems and, thus, to some extent, alleviates the curse of dimensionality. Second, incorporating additional terms into the learning problem penalizing the violation of the dynamic programming principles~\eqref{eq:dynamicalprog},\linebreak allows to compute a good approximation~$V^\eps_\theta$ of the optimal value function on the fly. As stated initially, the present example should be understood as a proof of concept and, following these first promising results, we believe that this approach to feedback learning deserves further investigations, both, from the theoretical and the numerical side. For example, it would be interesting to explore systematic ways of choosing the penalty parameters~$\gamma_1,\gamma_2$. However, this goes beyond the scope of the current paper and is left for future work.
\begin{table}[!htb]
\caption{Results on training set i.e.~$Y^{ad}_0=\bfy^t_0$.}\label{tab:val1}
\centering
\begin{tabular}{l c c c c c}

\toprule $\qquad \text{Penalty}$ & $\Err_{\clJ}$ & $\Err_{Y}$ & $\Err_{P}$ & $\Err_{U}$\\
\midrule
$\gamma_1= 0, \gamma_2=0$ & $ 0.15 \%$ & $0.04 \%$& $ 0.12 \%$ & $ 2.4 \%$ \\
\midrule
$\gamma_1= 0.1, \gamma_2=0.1$ & $0.36 \%$ & $0.1 \%$ & $0.24 \%$ & $5.5 \%$
\\
\midrule
$\gamma_1= 0.1, \gamma_2=0$ & $0.29 \%$ & $0.1 \%$ & $0.85 \%$ & $ 4.4 \%$
\\
\midrule
$\gamma_1= 1, \gamma_2=1$ & $ 0.64 \%$ & $ 0.25 \%$& $1 \%$ & $8.65 \%$ \\
\midrule
$\gamma_1= 0, \gamma_2=1$ & $ 0.1 \%$ & $ 0.05 \%$& $0.26 \%$ & $2.1 \%$ \\
\bottomrule
\end{tabular}
\\[0.5em]
\begin{tabular}{l c c c c c}

\toprule $\qquad \text{Penalty}$ & $\Err_{J}$ & $\Err_{V}$ & $\Err_{\partial V}$ & $d(V^\eps_\theta;V^*)$ & $d(\partial_y V^\eps_\theta;\partial_y V^*)$ \\
\midrule
$\gamma_1= 0, \gamma_2=0$ & $ 0.0003 \%$ & $ 79 \%$& $33 \%$ & $78.8 \%$ & $33.5 \%$ \\
\midrule
$\gamma_1= 0.1, \gamma_2=0.1$ & $ 0.001 \%$ & $0.03 \%$ & $7.4 \%$ & $0.03 \%$ & $7 \%$
\\
\midrule
$\gamma_1= 0.1, \gamma_2=0$ & $ 0.001 \%$ & $0.02 \%$ & $12.5 \%$ & $0.02 \%$ & $12.1 \%$
\\ \midrule
$\gamma_1= 1, \gamma_2=1$ & $ 0.005 \%$ & $ 0.007 \%$& $4.5 \%$ & $0.01 \%$ & $3.5 \%$
\\ \midrule
$\gamma_1= 0, \gamma_2=1$ & $ 0.003 \%$ & $ 88.8 \%$& $6.4 \%$ & $88.5 \%$ & $6.4 \%$\\
\bottomrule
\end{tabular}
\end{table}

\begin{table}[!htb]
\caption{Results on validation set i.e.~$Y^{ad}_0=\bfy^v_0$.}
\centering
\medskip
\begin{tabular}{l c c c c c}

\toprule $\qquad \text{Penalty}$ & $\Err_{\clJ}$ & $\Err_{Y}$ & $\Err_{P}$ & $\Err_{U}$ \\
\midrule
$\gamma_1= 0, \gamma_2=0$ & $ 0.23 \%$ & $0.06 \%$& $ 0.57 \%$ & $ 4.42 \%$ \\
\midrule
$\gamma_1= 0.1, \gamma_2=0.1$ & $0.51 \%$ & $0.15 \%$ & $1.1 \%$ & $9.2 \%$
\\ \midrule
$\gamma_1= 0.1, \gamma_2=0$ & $0.47 \%$ & $0.16 \%$ & $2.8 \%$ & $8.7 \%$
\\ \midrule
$\gamma_1= 1, \gamma_2=1$ & $ 0.85 \%$ & $ 0.35 \%$ & $ 6\%$ & $13.5 \%$
\\ \midrule
$\gamma_1= 0, \gamma_2=1$ & $ 0.25 \%$ & $ 0.1 \%$ & $ 1.3\%$ & $4.8 \%$\\
\bottomrule
\end{tabular}
\\[0.5em]
\begin{tabular}{l c c c c c}

\toprule $\qquad \text{Penalty}$ & $\Err_{J}$ & $\Err_{V}$ & $\Err_{\partial V}$ & $d(V^\eps_\theta;V^*)$ & $d(\partial_y V^\eps_\theta;\partial_y V^*)$ \\
\midrule
$\gamma_1= 0, \gamma_2=0$ & $ 0.002 \%$ & $ 78.7 \%$& $33.6 \%$ & $78.6 \%$ & $33.9 \%$ \\
\midrule
$\gamma_1= 0.1, \gamma_2=0.1$ & $0.007 \%$ & $0.03 \%$ & $8.9 \%$ & $0.03 \%$ & $7.9 \%$ \\
\midrule
$\gamma_1= 0.1, \gamma_2=0$ & $0.008 \%$ & $0.02 \%$ & $15.1 \%$ & $0.02 \%$ & $13.1 \%$ \\ \midrule
$\gamma_1= 1, \gamma_2=1$ & $ 0.02 \%$ & $0.009 \%$ & $ 9.8 \%$ & $ 0.01 \%$ & $ 4.2 \%$ \\
\midrule
$\gamma_1= 0, \gamma_2=1$ & $ 0.002 \%$ & $ 88 \%$ & $ 8.6 \%$ & $ 88 \%$ & $ 7.1 \%$ \\
\bottomrule
\end{tabular}
\label{tab:val2}
\end{table}

\appendix\label{app}
\def\theequation{\Alph{section}\arabic{equation}}
\setcounter{equation}{0}
\def\thetable{\Alph{section}\arabic{table}}

\section{Condition~\eqref{eq:kk10}} \label{app1}

Here we address condition~\eqref{eq:kk10}. Define $N_\epsilon=\lceil \frac{2\widehat M}{\eps} \rceil,\,\tilde M:=\eps N_\eps$ and introduce the equidistant grid $G=\{-\tilde M, (1-N_\eps)\eps, \,\dots,\,-\eps,0,\eps,\,\dots,\,(N_\eps-1)\eps, \tilde M\}$. Next endow the hypercube $[-\tilde M, \tilde M]^{n+1}$ with the $(n+1)-$ dimensional product of the grid $G$. These grid points define $\{Q_i\}_{i=1}^{(2N_\eps)^{n+1}}$ closed subhypercubes of dimension $\eps^n$ whose union covers $\bar K=[0,T] \times \bar B_{2\tilde M}(0)$.

We extend this $n+1$-dimensional grid by adding $k\ge \lceil\frac{1}{2}\sqrt{n}\rceil+1$ layers (again all of dimension~$\eps^n)$, to the surfaces of the preexisting grid, resulting in $\tilde N_\eps=(2N_\eps +2k)^{n+1}$ hypercubes whose union covers $ [- \tilde M-k\eps,\tilde M+k\eps]^{n+1}$. The subhypercubes are ordered in such a manner that the interiors ones $\{Q_i\}_{i=1}^{(2N_\eps +2(k-1))^{n+1}}$ are assembled first and the ones with a boundary face
\[
\mathop{\{ Q_i \}}\limits_{i=(2N_\eps +2(k-1)+1)^{n+1}}^{(2N_\eps +2k)^{n+1}}
\]
come last. The set of indices corresponding to interior hypercubes are denoted by $\clI$, those to boundary hypercubes by $\F$.

Next we introduce a staggered grid and place a node $x_i=(t_i,y_i)$ at the barycenter of each of the $Q_i,\, i=1,\,\dots,\,(2N_\eps +2k)^{n+1}$. We shall use the standard mollifier of radius $r_\eps$ defined by
\begin{equation*}
\psi(x)= \left\{
\begin{array}{ll}
\exp\left(\frac{1}{|\frac{x}{r_\eps}|^2-1}\right), & \text{ for} |x|\le r_\eps \\[1.4ex] 0, & \text{ for} |x|\le r_\eps,
\end{array}
\right.
\end{equation*}
where $r_\eps=\eps(\frac{1}{2}\sqrt{n}+.1)$. Note that by adding $.1$ in the previous expression the cube $[-\frac{\eps}{2},\frac{\eps}{2}]^n$ is contained in the interior of the support of $psi$. Finally we introduce $\psi_j(x)= \psi(x-x_j)$, for ${j\in \clI\cup \F}$ and
\begin{equation*}
\varphi_j= \frac{\psi_j}{\sum_{i\,\in\,\clI\,\cup\,\F} \psi_i}, \text{ for } j \in \clI.
\end{equation*}
Let us deduce the following properties:
\begin{enumerate}\romanenumi
\item\label{AppA.i} For each ${j\in \clI\cup \F}$ we have $supp\, \psi_j=\bar K_j$ where $K_j = \{x: |x-x_j|< r_\eps\}$.
\item\label{AppA.ii} By construction there exists $\frkm $ such that $\text{card} \{j: \psi_j(x) \neq 0 \} \le \frkm, \quad \forall x \in \R^{n+1}, \text{ and } \forall {j\in \clI\cup \F}$, for each $\eps\in (0,\eps_0]$.
\item\label{AppA.iii} For each $j\in \clI$ the denominator in the definition of $\varphi_j$ is different from zero. Hence $\varphi_j$ is well-defined with $supp\, \varphi_j=\bar K_j$ for $j\in \clI$ and $\varphi_j \colon\R^{n+1} \to [0,1]$, and it is $\clC ^\infty$ smooth.
\item\label{AppA.iv} $\bar{K} \subset \bigcup_{j\,\in\,\clI} Q_j\subset \bigcup_{j\,\in\,\clI} K_j$.
\item\label{AppA.v} $\text{card} \{j: \varphi_j(x) \neq 0 \} \le \frkm, \quad \forall x \in \R^{n+1}, \text{ and } \forall {j\in \clI}$, for each $\eps\in (0,\eps_0]$.
\item\label{AppA.vi} Due to the choice of $r_\eps$ the functions $\psi$ are uniformly bounded from below on $Q_j$ for each $j$, independent of $\eps\in (0,\eps_0]$. Moreover due to the boundedness of $\psi$ and by the definition of $\frkm$, there exists $\nu>0$ such that
\begin{equation*}
\sum_{i\,\in\,\clI\,\cup\,\F} \psi_i(x) \ge \nu, \forall x \in Q_j, \text{ with } j \in \clI\cup \F,
\end{equation*}
and thus in particular $\sum_{i\,\in\,\clI\,\cup\,\F} \psi_i(x) \ge \nu, \forall x \in \bar K$.
\item\label{AppA.vii}
\[
\supp\; \varphi_j \cap \bar K = \emptyset\; \forall j \in \F.
\]
This is a consequence of the fact that for $j \in \F$ we have $dist(x_j, \partial ([-\tilde M, \tilde M]^n))= \eps[(k-1)+\frac{1}{2}]$ and thus $dist(\partial K_j, \partial ([-\tilde M, \tilde M]^n)) \le \eps[(k-1)+\frac{1}{2}-r_\eps] = k-\frac{1}{2}(1+\sqrt{n})-.1>\eps (k-1 - \frac{1}{2}\sqrt{n})\linebreak>0$.
\item\label{AppA.viii} $\sum_{i\,\in\,\clI}\varphi_i =1, \; \forall x \in \bar K$. This is a consequence of~\eqref{AppA.vii} and the definition of $\varphi_j$.
\item\label{AppA.ix} $\|D^{j}\varphi_i\|_{C(\bar K_i\cap \bar K)} \le\bar \mu \eps^{-j}, \text{ for some } \bar \mu \text{ independent of } i\in \clI, \; \text{ and } j\in\{1,2\}.$
\end{enumerate}
Once we have verified $(ix)$, all the properties demanded in~\eqref{eq:kk10} on the partition of unity $\{\varphi_i\}_{i\in \clI}$ subordinate to~$K_i$ will be satisfied.

In the following calculations we repeatedly use that $ \nabla \sum _{i\in \clI}\varphi_i(x)=0$ for $x\in \bar K$. This follows from~\eqref{AppA.viii}. As short calculation shows that for each $j\in \clI$, each $x\in \bar K$, and $k,\ell \in \{1,\,\dots,\,n\}$
\begin{equation*}
\partial_{x_k} \varphi_j(x) = \frac{\partial_{x_k} \psi_j(x)}{\sum_{i\,\in\,\clI} \psi_i(x)}, \quad \partial_{x_\ell} \partial_{x_k} \varphi_j(x)= \frac{\partial_{x_\ell} \partial_{x_k} \varphi_j(x) \sum_{i\,\in\,\clI} \psi_i(x) - \psi_j \sum_{i\,\in\, \clI} \partial_{x_\ell} \partial_{x_k} \psi_i(x) }{\left(\sum_{i\,\in\,\clI} \psi_i(x)\right)^2 },
\end{equation*}
where we use that $\partial_{x_k}\sum_{i\,\in\,\clI} \psi_i(x)=0$ for $x\in \bar K$.

To obtain the required estimates we introduce for $\eta>0$
\begin{equation*}
\psi_\eta(x)= \left\{
\begin{array}{ll}
\exp\left(\frac{1}{\left|\frac{x}{\eta}\right|^2-1}\right), & \text{ for} |x|\le \eta \\[1.4ex] 0, & \text{ for} |x|\le \eta.
\end{array}
\right.
\end{equation*}
Then we have
\begin{align*}
\partial_{x_k} \psi_\eta &=- \psi_\eta \, \frac{2x_k}{\eta^2 \left(\left|\frac{x}{\eta}\right|^2-1\right)^2},
\\
\left(\partial_{x_k}\right)^2 \psi_\eta &=
\frac{2\psi_\eta}{\eta^2\left(\left|\frac{x}{\eta}\right|^2-1\right)^4}\left[\frac{2x_k^2}{\eta^2}-\left(\left|\frac{x}{\eta}\right|^2-1\right)^2 +\frac{4x_k^2}{\eta^2}\left(\left|\frac{x}{\eta}\right|^2-1\right)\right],
\end{align*}
and for $k\ne \ell$
\begin{equation*}
\partial_{x_\ell} \partial_{x_k}\psi_\eta=
\frac{2\psi_\eta\,x_\ell x_k}{\eta^2\left(\left|\frac{x}{\eta}\right|^2-1\right)^4} \left[\frac{2}{\eta^2} + \frac{4}{\eta^2}\left(\left|\frac{x}{\eta}\right|^2 -1\right)\right]=
\frac{2\psi_\eta\, x_\ell x_k}{\eta^2\left(\left|\frac{x}{\eta}\right|^2-1\right)^4} \left[\frac{-2}{\eta^2} + \frac{4}{\eta^4}|x|^2 \right].
\end{equation*}
Considering the behavior of $\partial_{x_k} \psi_\eta$ and $\partial_{x_\ell} \partial_{x_k}\psi_\eta$ separately on the ball $B_{\frac{\eta}{2}}(0)$ and its complement in $B_{{\eta}}(0)$, it follows that these functions behave like $O(\frac{1}{\eta})$ and $O(\frac{1}{\eta^2})$. Applying these estimates in the expressions for the first and second derivatives for $\varphi_j$ and using the lower bound established in $(vi)$ we obtain $(ix)$.



\section{Perturbation results} \label{app:perturbation}
Here we collect pertinent existence and stability results for dynamical systems. The constant $\widehat{M}(0)$ appearing below relates to Assumption $(\textbf{A.2})$.

\begin{prop}\label{prop:s}
Let~$\bfy \in \clC (Y_0;W_T)$,~$\bfy(y_0) \in \clY _{ad}$ for all~$y_0 \in Y_0$,~$\delta \bfv  \in \clC (Y_0;L^2(I;\R^n))$, as well as~$\delta \bfy_0 \in \clC (Y_0;\R^n)$ be given. Moreover let~$A \colon I \times \R^n \to \R^{n \times n} $ be continuous, and denote by~$\bfA  \colon L^\infty(I;\R^n) \to \clB (L^2(I;\R^n))$ the induced Nemitsky operator i.e.
\begin{align*}
\bfA (y) \delta y = A(t,y(t)) \delta y(t) \quad \forall \delta y \in L^2\left(I;\R^n\right),~y \in L^\infty\left(I;\R^n\right)
\end{align*}
and a.e.~$t \in I$. Then there is~$\delta \bfy \in \clC (Y_0;W_T)$ such that~
\begin{align}\label{eq:linearensemble}
\dot{\delta y}= \bfA (y)\delta y+v,~\delta y (0)= \delta y_0
\end{align}
for~$y\coloneqq \bfy(y_0)$,~$\delta v\coloneqq \delta \bfv (y_0) $,~$\delta y_0 \coloneqq \delta \bfy(y_0)$ and all~$y_0 \in Y_0$. It satisfies
\begin{align}\label{eq:aprioriindependent}
\left\|\delta \bfy(y_0)\right\|_{W_T} \leq C \left(\left\|\delta \bfv (y_0)\right\|_{L^2}+\left|\delta \bfy_0(y_0)\right| \right)
\end{align}
for some~$C>0$ depending continuously on $\max_{(\tau,y)\in I \times \bar{B}_{2 \widehat{M}(0)} }\|A(\tau,y)\|_{\R^{n \times n}}$, and independent of $y_0\in Y_0$.
\end{prop}

\begin{proof}
Let~$y_0 \in Y_0$ be arbitrary but fixed. Then there is a unique solution~$\delta y \in W_T$ to~\eqref{eq:linearensemble} which satisfies
\begin{align*}
\frac{1}{2} \left|\delta y (t)\right|^2 &= \frac{1}{2} \left|\delta y (0)\right|^2+ \int^t_0 \left(\dot{\delta y}(s),\delta (y)(s)\right)~\de s \\
&= \frac{1}{2} |\delta y_0|^2+ \int^t_0 \left(\delta y(s), A(t,y(t))\delta (y)(s)\right)+ \left(\delta v(s), \delta y(s)\right)~\de s \\
& \leq \frac{1}{2} |\delta y_0|^2+ \frac{1}{2} |\delta v|^2_{L^2}+ \frac{1}{2}\int^t_0 \left(2 \max_{(\tau,y)\in I \times \bar{B}_{2 \widehat{M}(0)} }\|A(\tau,y)\|_{\R^{n \times n}} +1\right) |\delta y(s)|^2~\de s
\end{align*}
for all~$t \in I$. Setting
\begin{align*}
L \coloneqq \left(2 \max_{(z,y)\in I \times \bar{B}_{2 \widehat{M}(0)} }\|A(z,y)\|_{\R^{n \times n}} +1\right),
\end{align*}
Gronwall's inequality implies that
\begin{align*}
\|\delta y\|_{L^\infty} \leq e^{TL}\left(|\delta y_0|+ |\delta v|_{L^2}\right).
\end{align*}
By~\eqref{eq:linearensemble} we further get $
\|\dot{\delta y}\|_{L^2} \leq L (\|\delta y\|_{L^2} + \|v\|_{L^2}), $ which implies~\eqref{eq:aprioriindependent}. Next, let~$y^k_0 \in Y_0$ denote a convergent sequence with limit~$y_0$. For abbreviation set
\begin{align*}
\delta y_k \coloneqq \delta \bfy\left(y^k_0\right),~y_k\coloneqq \bfy\left(y^k_0\right),~\delta v_k \coloneqq \delta \bfv \left(y^k_0\right),~\delta y_0 \coloneqq \delta \bfy_0\left(y^k_0\right)
\end{align*}
as well as
\begin{align*}
y\coloneqq \bfy(y_0),~\delta v \coloneqq \delta \bfv (y_0),~\delta y_0 \coloneqq \delta \bfy_0(y_0).
\end{align*}
Note that~$\delta y_k$ is uniformly bounded in~$W_T$ by~\eqref{eq:aprioriindependent}. Thus it admits a subsequence, denoted by the same index, with~$\delta y_k \rightharpoonup \delta y$ in~$W_T$ for some~$\delta y \in W_T$. This implies
\begin{align*}
\delta y_k (0) \rightarrow \delta y(0)~\text{in}~\R^n,~\delta y_k \rightarrow \delta y~\text{in}~L^\infty\left(I;\R^n\right),~\dot{\delta y}_k \rightharpoonup \dot{\delta y}~\text{in}~L^2\left(I;\R^n\right).
\end{align*}
Moreover, due to the continuity of~$\bfy,\delta \bfv $ and~$\delta \bfy_0$, we get
\begin{align*}
\dot{\delta y}_k=\bfA (y_k)\delta y_k+ \delta v_k \rightarrow \bfA (y)\delta y+ \delta v~\text{in}~L^2\left(I;\R^n\right),~\delta y_k (0) \rightarrow \delta y_0~\text{in}~\R^n.
\end{align*}
Summarizing the previous observations we conclude that
\begin{align*}
\dot{\delta y}= \bfA (y)\delta y+\delta v,~\delta y(0)= \delta y_0
\end{align*}
as well as~$\delta y_k \rightarrow \delta y$ in~$W_T$, and thus~$\delta y= \delta \bfy(y_0) $. By uniqueness of solutions to the above equation~$\delta \bfy(y^k_0)\rightarrow \delta \bfy(y_0) $ for the whole sequence in~$W_T$ follows, and therefore~$\delta \bfy \in \clC (Y_0;W_T)$.
\end{proof}

Next we address nonlinear systems of the form:
\begin{align}\label{eq:pertstate}
\dot{y_v}= \bff(y_v)+ \bfg(y_v) \F^*(y_v)+v, \quad y_v(0)=y_0
\end{align}
where~$v \in L^2(I ; \R^n)$ is a perturbation.

\begin{prop}\label{thm:existspert}
Let Assumption~\ref{ass:feedbacklaw} hold. Then there exist an open neighbourhood~$V_1 \subset L^2(I; \R^n)$ of~$0$ and an open neighbourhood~$\bfy_0$ of $Y_0$ such that~\eqref{eq:pertstate} admits a unique solution~$y_v=\bfy^v (y_0)\in \clY _{ad} $ for every pair~$(v,y_0)\in V_1 \times
\bfy_0$. Moreover the mapping
\begin{align}
\bfy^\bullet(\bcdot) \colon V_1 \times \bfy_0 \to
\clY _{ad}, \quad (v, y_0) \mapsto \bfy^v(y_0)
\end{align}
is continuously Fr\'echet differentiable.
\end{prop}
\begin{proof}
Define the mapping
\begin{align*}
G \colon \clY _{ad} \times \R^n \times L^2\left(I; \R^n\right)
\to L^2\left(I; \R^n\right) \times \R^n
\end{align*}
with
\begin{align*}
G(y,y_0,v)= \left(
\begin{array}{c}
\dot{y}-\bff(y)-\bfg(y)\F^*(y)-v\\
y(0)-y_0\\
\end{array}
\right).
\end{align*}
Now fix an arbitrary~$\bar{y}_0 \in Y_0$ and, utilizing $(\mathbf{A.3})$ denote by $\bar{y}=\bfy^*(y_0)\in \intt
\clY _{ad}$ the unique solution in~$\clY _{ad}$ to the unperturbed closed loop system~$G(\bar{y},\bar{y}_0,0)=0$. Since~$G$ is of class~$\clC ^1$ in a neighborhood of~$(\bar{y},\bar{y}_0,0)$ we have
\begin{align*}
D_y G(y,y_0,v)\delta y = \left(
\begin{array}{c}
\dot{\delta y}-D\bff(y)\delta y-\left\lbrack D\bfg(y) \delta y \right\rbrack
\F^*(y) - \bfg(y) \partial_y \F^*(y)\delta y\\
\delta y(0)\\
\end{array}
\right).
\end{align*}
It is straightforward that the linearized equation
\begin{align*}
D_y G(\bar{y},\bar{y}_0,v)\delta y=\left(
\begin{array}{c}
\delta v \\ \delta y_0\\
\end{array}
\right)
\end{align*}
admits a unique solution~$\delta \bar{y}\in W_T$ for every~$\delta v \in L^2(I;\R^n),~\delta y_0 \in \R^n$. Moreover, applying Gronwall's lemma yields~$c>0$ independent of~$\bar{y},~\bar{y}_0$ with
\begin{align*}
\wnorm{\delta \bar{y}} \leq c(\|\delta v\|_{L^2(I;\R^n)}+|\delta y_0|), \quad \forall \delta v \in L^2\left(I;\R^n\right),~\delta y_0 \in \R^n.
\end{align*}
Thus from the implicit function theorem we get constants~$\kappa_1=\kappa_1(\bar{y}_0)$ and~$\kappa_2=\kappa_2(\bar{y}_0)$, such that for every~$y_0\in \R^n$ with~$|y_0-\bar{y}_0|<\kappa_1$ and~$\|v\|_{L^2(I;\R^n)}< \kappa_2$ there exists~$\bfy^v(y_0) \in \clY _{ad}$ with~$G(\bfy^v(y_0),y_0,v)=0$. By $(\textbf{A.1})$ it is the unique solution to~\eqref{eq:pertstate} in~$\clY _{ad}$. Moreover, the mapping
\begin{align*}
\bfy^\bullet(\cdot) \colon B_{\kappa_2}(0) \times B_{\kappa_1}(\bar{y}_0) \to \clY , \quad (v,y_0) \mapsto
\bfy^v(y_0)
\end{align*}
is of class~$\clC ^1$. Observe that repeating this argument for every~$y_0 \in Y_0$ yields an open covering of~$Y_0$ i.e.
\begin{align*}
Y_0 \subset \bigcup_{\bar{y}_0\,\in\,Y_0} B_{\kappa_1(\bar{y}_0)}(\bar{y}_0).
\end{align*}
Since~${Y}_0$ is compact there exists a finite set of initial conditions~$\{\bar{y}^i_0\}^N_{i=1}\subset Y_0$, including~$0$, such that
\begin{align*}
Y_0 \subset \bfy_0 :=\bigcup^N_{i=1} B_{\kappa_1\left(\bar{y}^i_0\right)}\left(\bar{y}^i_0\right).
\end{align*}
Set~$V= \bigcap^N_{i=1} B_{\kappa_2(\bar{y}^i_0)}(0) \subset L^2(I;\R^n) $. Summarizing these arguments yields the existence of a~$\clC ^1$-mapping
\begin{align*}
\bfy^\cdot(\cdot) \colon V \times \bfy_0 \to
\clY _{ad}, \quad \bfy^v(y_0)~\text{ uniquely solves}~\eqref{eq:pertstate}~\text{in}~\clY_{ad}.
\end{align*}
\end{proof}

We use the following consequences of the previous proposition.

\begin{enonce}
{Corollary} \label{coroll:locallipschitzofstate}
There exists an open neighborhood~$V_2 \subset V_1 \subset L^2(I;\R^n)$ of~$0$ as well as~$c>0$ such that
\begin{align*}
\left\|\bfy^{v_1}(y_0)-\bfy^{v_2}(y_0)\right\|_{W_T}
 \leq c
\|v_1-v_2\|_{L^2(I;\R^n)} \quad \forall y_0 \in Y_0,~v_1 \in V_2, v_2 \in V_2
\\
\intertext{and}
\left\|\bfy^v(y_0)\right\|_{W_T} \leq M_{Y_0} +c \|v\|_{L^2(I;\R^n)} \quad
\forall y_0 \in Y_0,~v \in V_2,
\end{align*}
hold. Here~$M_{Y_0}$ denotes the constant from~\eqref{A.3}.
\end{enonce}

\begin{proof}
The first assertion follows from the continuous differentiability of $v\to \bfy^v(y_0)$ and compactness of $Y_0$. To verify the second we use that $\bfy^*(y_0) = \bfy^0(y_0)$ and estimate
\begin{align*}
\left\|\bfy^v(y_0)\right\|_{W_T} &\leq \left\|\bfy^*(y_0)\right\|_{W_T} +
\left\|\bfy^v(y_0)-\bfy^0(y_0)\right\|_{W_T}.
\end{align*}
The claim now follows from the first inequality and~\eqref{A.3}.
\end{proof}

\section*{Conflicts of interest}
The authors have no conflict of interest to declare.

\bibliographystyle{crunsrt}
\bibliography{crmeca20221181}
\end{document}
