remarkRemark exampleExample hypothesisHypothesis assumptionAssumption
claimClaim ex_supplement An Optimal Control Problem with Complementarity
Constraints J. LUO and X. CHEN

# An Optimal Control Problem with Terminal Stochastic Linear Complementarity
Constraints††thanks: Submitted to the editors DATE. This work is supported by
the Hong Kong Research Grants Council grant PolyU15300021 and the CAS AMSS-
PolyU Joint Laboratory in Applied Mathematics.

Jianfeng Luo CAS AMSS-PolyU Joint Laboratory of Applied Mathematics
(Shenzhen), The Hong Kong Polytechnic University Shenzhen Research Institute,
Shenzhen, China<EMAIL_ADDRESS>Xiaojun Chen Department of Applied
Mathematics, The Hong Kong Polytechnic University, Hong Kong, China ().
<EMAIL_ADDRESS>
###### Abstract

In this paper, we investigate an optimal control problem with terminal
stochastic linear complementarity constraints (SLCC), and its discrete
approximation using the relaxation, the sample average approximation (SAA) and
the implicit Euler time-stepping scheme. We show the existence of feasible
solutions and optimal solutions to the optimal control problem and its
discrete approximation under the conditions that the expectation of the
stochastic matrix in the SLCC is a Z-matrix or an adequate matrix. Moreover,
we prove that the solution sequence generated by the discrete approximation
converges to a solution of the original optimal control problem with
probability 1 as $\epsilon\downarrow 0$, $\nu\to\infty$ and $h\downarrow 0$,
where $\epsilon$ is the relaxation parameter, $\nu$ is the sample size and $h$
is the mesh size. We also provide asymptotics of the SAA optimal value and
error bounds of the time-stepping method. A numerical example is used to
illustrate the existence of optimal solutions, the discretization scheme and
error estimation.

###### keywords:

linear complementarity problem, sample average approximation, implicit Euler
time-stepping, convergence analysis, error estimate.

49M25, 49N10, 90C15, 90C33

##  1 Introduction

Let $\xi$ be a random variable defined in the probability space
$(\Omega,\mathcal{F},\mathcal{P})$ with support set
$\Xi:=\xi(\Omega)\subseteq\mathbb{R}^{d}$. Let $\|\cdot\|$ denote the
Euclidean norm of a vector and a matrix. Let $L^{2}(0,T)^{n}$ denote the
Banach space of all quadratically Lebesgue integrable functions mapping from
$(0,T)$ to $\mathbb{R}^{n}$, which is equipped with the norm

$\|x\|_{L^{2}}:=\left(\int_{0}^{T}\|x(t)\|^{2}dt\right)^{\frac{1}{2}},\,\,\,\forall\,\,x\in
L^{2}(0,T)^{n}.$

Denote $H^{1}(0,T)^{n}$ the space of all functions $x\in L^{2}(0,T)^{n}$ whose
components $x_{1},\cdot\cdot\cdot,x_{n}:(0,T)\rightarrow\mathbb{R}$ possess
weak derivatives such that the function $\dot{x}\in L^{2}(0,T)^{n}$. A
suitable norm in $H^{1}(0,T)^{n}$ is defined by

$\|x\|_{H^{1}}:=\left(\|x\|_{L^{2}}^{2}+\|\dot{x}\|_{L^{2}}^{2}\right)^{\frac{1}{2}},\,\,\,\forall\,\,x\in
H^{1}(0,T)^{n}.$

In this paper, we aim to find an optimal solution $(x,u)\in
H^{1}(0,T)^{n}\times L^{2}(0,T)^{m}$ of the following optimal control problem
with terminal stochastic linear complementarity constraints:

(1)
$\displaystyle\min_{x,u}\,\,\Phi(x,u):=\mathbb{E}[F(x(T),\xi)]+\frac{1}{2}\|x-x_{d}\|_{L^{2}}^{2}+\frac{\delta}{2}\|u-u_{d}\|_{L^{2}}^{2}$
$\displaystyle\textrm{s.t.}\left\\{\begin{aligned} &\left.\begin{aligned}
&\dot{x}(t)=Ax(t)+Bu(t),&\\\ &Cx(t)+Du(t)-f(t)\leq 0,&\\\
\end{aligned}\right\\}a.e.\,\,t\in(0,T),&\\\ &0\leq
x(T)\perp\mathbb{E}[M(\xi)x(T)+q(\xi)]\geq 0,&\\\
&x(0)=x_{0},\,\mathbb{E}[g(x(T),\xi)]\in K,\end{aligned}\right.$

where $A\in\mathbb{R}^{n\times n}$, $B\in\mathbb{R}^{n\times m}$,
$C\in\mathbb{R}^{l\times n}$, $D\in\mathbb{R}^{l\times m}$,
$x_{0}\in\mathbb{R}^{n}$ and $f\in L^{2}(0,T)^{l}$ are given, $\delta>0$ is a
scalar, $K\subseteq\mathbb{R}^{\ell}$ is a nonempty, closed and convex set.
Here $x_{d}\in L^{2}(0,T)^{n}$ and $u_{d}\in L^{2}(0,T)^{m}$ are the given
desired state and control, respectively. The mappings
$F:\mathbb{R}^{n}\times\Xi\rightarrow\mathbb{R}$ and
$g:\mathbb{R}^{n}\times\Xi\rightarrow\mathbb{R}^{\ell}$ are continuously
differentiable over $\mathbb{R}^{n}$ for almost every $\xi\in\Xi$, and are
measurable on $\Xi$ for each $z\in\mathbb{R}^{n}$. Two mappings
$M:\Xi\rightarrow\mathbb{R}^{n\times n}$ and $q:\Xi\rightarrow\mathbb{R}^{n}$
are also measurable in $\xi\in\Xi$.

In [2], Benita and Mehlita studied an optimal control problem with terminal
deterministic complementarity constraints, which has many interesting
practical applications in multi-agent control networks. They derived some
stationarity conditions and presented constraint qualifications which ensure
that these conditions hold at a local optimal solution of the optimal control
problem under the assumption that the feasible set is nonempty. However,
sufficient conditions were not given for the existence of $x(T)$ such that the
terminal deterministic complementarity constraints

(2) $\displaystyle 0\leq x(T)\perp\bar{M}x(T)+\bar{q}\geq 0,$ (3)
$\displaystyle\bar{g}(x(T))\in K,$

hold, where $\bar{M}\in\mathbb{R}^{n\times n}$, $\bar{q}\in\mathbb{R}^{n}$ and
$\bar{g}:\mathbb{R}^{n}\rightarrow\mathbb{R}^{\ell}$. Problem (1) extends the
problem in [2] to terminal stochastic case in uncertain environment.

Optimal control with differential equations and complementarity constraints
provides a powerful modeling paradigm for many practical problems such as the
optimal control of electrical networks with diodes and/or MOS transistors [4]
and dynamic optimization of chemical processes [22]. It can also be derived
from the KKT conditions of a bilevel optimal control if the lower level
problem is convex and satisfies a constraint qualification [17]. A series of
works [5, 6, 10, 14, 24] are devoted to the study of optimal control problems
with complementarity constraints. It should be noted that these papers focus
on deterministic problems, where the system coefficients including system
parameters and boundary/initial conditions are perfectly known. On the other
hand, optimal control problems with stochastic differential equation
constraints under uncertain environment have been extensively studied [16, 19,
20]. These papers investigate theory and algorithms for optimal control when
the parameters in the differential equations have noise and uncertainties.
However, there is very little research on optimal control with terminal
stochastic complementarity constraints.

The main contributions of this paper are summarized as follows. We show the
existence of feasible solutions to the optimal control problem (1) under the
conditions that $\mathbb{E}[M(\xi)]$ is a Z-matrix or an adequate matrix,
which gives reasonable conditions for the existence of $x(T)$ such that
(2)-(3) hold. Moreover, we prove the existence of feasible solutions and
optimal solutions to the discrete approximation using the relaxation, the
sample average approximation (SAA) and the implicit Euler time-stepping scheme
under the same conditions. In the convergence analysis, we prove that the
solution sequence generated by the discrete approximation converges to a
solution of the original optimal control problem with probability 1 (w.p.1) as
$\epsilon\downarrow 0$, $\nu\to\infty$ and $h\downarrow 0$, where $\epsilon$
is the relaxation parameter, $\nu$ is the sample size and $h$ is the mesh
size. We also provide asymptotics of the SAA optimal value and error bounds of
the time-stepping method. These results extend the approximation error of the
Euler time-stepping method of an optimal control problem with convex terminal
constraints to nonconvex terminal stochastic complementarity constraints.

The paper is organised as follows: Section 2 deals with the existence of
feasible solutions of problem (1). Section 3 studies the existence of feasible
solutions of the relaxation and the SAA of (1) and the convergence to the
original problem (1) as the relaxation parameter goes to zero and the sample
size approaches to infinity. In Section 4, we study the convergence of the
time-stepping scheme and show the convergence properties of the discrete
method using the SAA and the implicit Euler time-stepping scheme. A numerical
example is given in Section 5 to illustrate the theoretical results obtained
in this paper. Final conclusion remarks are presented in Section 6.

###  1.1 Notation and assumptions

Throughout this paper we use the following notation. For a matrix
$A\in\mathbb{R}^{m\times n}$, $A_{i\bullet}$ and $A_{\bullet j}$ denote its
$i$th row and $j$th column, respectively, $A^{\top}$ denotes its transpose
matrix, and $A^{\dagger}$ is its pseudoinverse matrix. If $A$ possesses full
row rank $m$, we have $A^{\dagger}=A^{\top}(AA^{\top})^{-1}$. Let $I$ denote
the identity matrix with a certain dimension.

For sets ${S}_{1},S_{2}\subseteq\mathbb{R}^{n}$, we denote the distance from
$v\in\mathbb{R}^{n}$ to $S_{1}$ and the deviation of the set $S_{1}$ from the
set $S_{2}$ by $\textrm{dist}(v,S_{1})=\inf_{v^{\prime}\in
S_{1}}\|v-v^{\prime}\|,$ and $\mathbb{D}(S_{1},S_{2})=\sup_{v\in
S_{1}}\textrm{dist}(v,S_{2})$, respectively. For sets ${S}_{1},S_{2}\subseteq
H^{1}(0,T)^{n}\times L^{2}(0,T)^{m}$, we denote the distance from
$(v_{1},v_{2})\in H^{1}(0,T)^{n}\times L^{2}(0,T)^{m}$ to $S_{1}$ by
$\textrm{dist}((v_{1},v_{2}),S_{1})=\inf_{(v^{\prime}_{1},v^{\prime}_{2})\in
S_{1}}(\|v_{1}-v^{\prime}_{1}\|_{H^{1}}+\|v_{2}-v^{\prime}_{2}\|_{L^{2}}),$
and the deviation of the set $S_{1}$ from the set $S_{2}$ by
$\mathbb{D}(S_{1},S_{2})=\sup_{(v_{1},v_{2})\in
S_{1}}\textrm{dist}((v_{1},v_{2}),S_{2})$. Let
$\mathcal{B}(v,\varepsilon)=\\{w:\|w-v\|\leq\varepsilon\\}$ be the closed ball
centered at $v$ with the radius of $\varepsilon$. Let int$S$ denote the
interior of a set $S$. Let $[N]=\\{1,2,\ldots,N\\}.$

We give two basic assumptions to ensure that problem (1) is well defined.

There exist four nonnegative measurable functions $\kappa_{i}(\xi)$ with
$\mathbb{E}[\kappa_{i}(\xi)]<\infty$ $(i=1,2,3,4)$ such that for any $z_{1}$,
$z_{2}\in\mathbb{R}^{n}$,

$|F(z_{1},\xi)-F(z_{2},\xi)|\leq\kappa_{1}(\xi)\|z_{1}-z_{2}\|,\,\,\|g(z_{1},\xi)\|\leq\kappa_{2}(\xi)\|z_{1}\|,\,\,a.e.\,\,\xi\in\Xi,$

and

$\|M(\xi)\|\leq\kappa_{3}(\xi)\,\,\,\textrm{and}\,\,\,\|q(\xi)\|\leq\kappa_{4}(\xi),\,\,\,\forall\,\,\xi\in\Xi.$

The matrix $D\in\mathbb{R}^{l\times m}$ is full row rank with $l<m$ and the
matrix

$\mathcal{R}:=[BY\,\,\,\,(A-BD^{\dagger}C)BY\,\,\,(A-BD^{\dagger}C)^{2}BY\,\,\,\cdot\cdot\cdot\,\,\,(A-BD^{\dagger}C)^{n-1}BY]\in\mathbb{R}^{n\times
n(m-l)}$

is also full row rank, where $Y\in\mathbb{R}^{m\times(m-l)}$ is a matrix with
full column rank $m-l$ such that $DY=0$.

##  2 Existence of optimal solutions of problem (1)

In this section, we first investigate the feasibility of problem (1). We call
$(x,u)\in H^{1}(0,T)^{n}\times L^{2}(0,T)^{m}$ a feasible solution of (1) if
it satisfies the constraints in (1).

For an index set $J\subseteq\\{1,2,\cdot\cdot\cdot,n\\}$, let $|J|$ denote its
cardinality and $J^{c}$ denote its complementarity set. We denote by
$q_{J}\in\mathbb{R}^{|J|}$ the subvector formed from a vector
$q\in\mathbb{R}^{n}$ by picking the entries indexed by $J$ and denote by
$M_{J_{1},J_{2}}\in\mathbb{R}^{|J_{1}|\times|J_{2}|}$ the submatrix formed
from a matrix $M\in\mathbb{R}^{n\times n}$ by picking the rows indexed by
$J_{1}$ and columns indexed by $J_{2}$. Let ${\cal
J}=\\{J\,:\,\mathbb{E}[M_{J,J}(\xi)]$ is nonsingular$\\}$ and

$\beta=\left\\{\begin{array}[]{ll}1&{\rm if}\,{\cal J}=\emptyset,\\\
\max\\{\|(\mathbb{E}[M_{J,J}(\xi)])^{-1}\|_{1}\,\,|\,\,J\in{\cal J}\\}&{\rm
otherwise.}\end{array}\right.$

A square matrix is said to be a P-matrix if all its principal minors are
positive. A square matrix is said to be a Z-matrix if its off-diagonal entries
are non-positive. A matrix $\mathbb{E}[M(\xi)]\in\mathbb{R}^{n\times n}$ is
called column adequate if for each $z\in\mathbb{R}^{n}$,
$z_{i}(\mathbb{E}[M(\xi)]z)_{i}\leq 0$ for all $i=1,2,\cdot\cdot\cdot,n$
implies $\mathbb{E}[M(\xi)]z=0$. The matrix $\mathbb{E}[M(\xi)]$ is row
adequate if $\mathbb{E}[M(\xi)]^{\top}$ is column adequate and it is adequate
if it is both column and row adequate [11]. It is known that a P-matrix is
adequate and a symmetric positive semi-definite matrix is also adequate [11,
Theorem 3.1.7, Theorem 3.4.4]. However, an adequate matrix may neither be a
P-matrix nor a positive semi-definite matrix [11].

######  Theorem 2.1.

Let Assumption 1.1 and Assumption 1.1 hold. Suppose the following three
conditions hold: (i)
$\mathcal{B}(0,\beta\mathbb{E}[\kappa_{2}(\xi)]\|\mathbb{E}[q(\xi)]\|_{1})\subseteq
K$, (ii) the set

$\mathcal{V}:=\\{v\in\mathbb{R}^{n}\,|\,\mathbb{E}[M(\xi)v+q(\xi)]\geq
0,\,v\geq 0\\}$

is nonempty, and (iii) $\mathbb{E}[M(\xi)]$ is an adequate matrix or a
Z-matrix. Then problem (1) has a feasible solution $(x,u)\in
H^{1}(0,T)^{n}\times L^{2}(0,T)^{m}$. Moreover, the optimal control problem
(1) admits an optimal solution if $\mathbb{E}[F(\cdot,\xi)]$ is bounded from
below.

######  Proof 2.2.

Following Lemma 7.2 in [2], Assumption 1.1 implies that the following system

(4) $\left\\{\begin{aligned} &\dot{x}(t)=Ax(t)+Bu(t),&\\\
&Cx(t)+Du(t)=p(t),&\\\ &x(0)=x_{0},\,\,\,x(T)=b,&\\\
\end{aligned}\right.\,\,a.e.\,\,t\in(0,T),$

admits a solution $(x,u)\in H^{1}(0,T)^{n}\times L^{2}(0,T)^{m}$ for any $p\in
L^{2}(0,T)^{l}$ and $b\in\mathbb{R}^{n}$. If we set $p(t)=f(t)+\tilde{p}(t)$
in (4) for arbitrary $\tilde{p}\in L^{2}(0,T)^{l}$ with $\tilde{p}(t)\leq 0$
and $f(t)$ in (1), then the following problem

(5) $\left\\{\begin{aligned} &\dot{x}(t)=Ax(t)+Bu(t),&\\\
&Cx(t)+Du(t)-f(t)\leq 0,&\\\ &x(0)=x_{0},\,\,\,x(T)=b,&\\\
\end{aligned}\right.\,\,a.e.\,\,t\in(0,T),$

has a solution $(x,u)\in H^{1}(0,T)^{n}\times L^{2}(0,T)^{m}$ for any
$b\in\mathbb{R}^{n}$.

Now we show the solution set of the following stochastic constrained LCP is
nonempty,

(6) $\left\\{\begin{aligned}
&\min\\{x(T),\,\mathbb{E}[M(\xi)x(T)+q(\xi)]\\}=0,&\\\
&\mathbb{E}[g(x(T),\xi)]\in K.\end{aligned}\right.$

Following Corollary 3.3.5 and Theorem 3.11.6 in [11], the LCP in (6) has a
solution from the assumption that the set $\mathcal{V}$ is nonempty and
$\mathbb{E}[M(\xi)]$ is adequate or a Z-matrix. Let $x^{*}(T)$ be a sparse
solution of the LCP in (6), which is defined as

$\displaystyle x^{*}(T)\in$ $\displaystyle\quad{\rm argmin}\,\|x(T)\|_{0}$
$\displaystyle\quad{\rm s.t.}\quad\quad 0\leq
x(T)\perp\mathbb{E}[M(\xi)x(T)+q(\xi)]\geq 0,$

where $\|x(T)\|_{0}=$ the number of nonzero components of $x(T)$. If there is
no $J$ such that $\mathbb{E}[M_{J,J}(\xi)]$ is nonsingular, that is, ${\cal
J}=\emptyset$, then by Theorem 2.2 in [8] and Theorem 3.1 in [9],
$\|x^{*}(T)\|_{0}=\|x^{*}(T)\|_{1}=0.$ Hence, we have

(7) $\|x^{*}(T)\|\leq\beta\|\mathbb{E}[q(\xi)]\|_{1}.$

If there is $J$ such that $x^{*}(T)_{J}>0$ and $x^{*}(T)_{J^{c}}=0$, from
Theorem 2.2 in [8] and Theorem 3.1 in [9], we know that
$\mathbb{E}[M_{J,J}(\xi)]$ is nonsingular and

$x^{*}(T)=-(I-\Lambda+\Lambda\mathbb{E}[M(\xi)])^{-1}\Lambda\mathbb{E}[q(\xi)],$

where $\Lambda$ is a diagonal matrix with $\Lambda_{i,i}=1,$ if $i\in J$ and
$\Lambda_{i,i}=0,$ if $i\in J^{c}$. Moreover, from

$\|(I-\Lambda+\Lambda\mathbb{E}[M(\xi)])^{-1}\Lambda\|\leq\max\\{\|(\mathbb{E}[M_{J,J}(\xi)])^{-1}\|_{1}\,|\,J\in{\cal
J}\\},$

we obtain (7) for ${\cal J}\neq\emptyset$.

Therefore, from Assumption 1.1 and assumption (i) of this theorem, we have

$\displaystyle\|\mathbb{E}[g(x^{*}(T),\xi)]\|\leq\mathbb{E}[\kappa_{2}(\xi)]\|x^{*}(T)\|\leq\mathbb{E}[\kappa_{2}(\xi)]\|x^{*}(T)\|_{1}\leq\beta\mathbb{E}[\kappa_{2}(\xi)]\|\mathbb{E}[q(\xi)]\|_{1},$

which implies that $\mathbb{E}[g(x^{*}(T),\xi)]\in K$. Hence the solution set
of (6) is nonempty.

Similar to the proof of Theorem 5.1 in [2], we can derive the existence of
optimal solutions to problem (1) under the assumption that
$\mathbb{E}[F(\cdot,\xi)]$ is bounded from below.

######  Remark 2.3.

The constrained LCP (6) may have multiple solutions or may not have a
solution. If $\mathbb{E}[M(\xi)]$ is a P-matrix, then for any
$\mathbb{E}[q(\xi)]$, the LCP in (6) has a unique solution $x(T)$. In such
case, if $\mathbb{E}[g(x(T),\xi)]\in K$, then (6) has a unique solution,
otherwise (6) does not have a solution. If $\mathbb{E}[M(\xi)]$ is a Z-matrix
or an adequate matrix, the LCP in (6) may have multiple solutions, while some
solutions can be bounded by $\beta\|\mathbb{E}[q(\xi)]\|_{1}$. When
$\mathcal{B}(0,\beta\mathbb{E}[\kappa_{2}(\xi)]\|\mathbb{E}[q(\xi)]\|_{1})$
$\subseteq K$, some solutions of the LCP satisfy $\mathbb{E}[g(x(T),\xi)]\in
K$ and thus the constrained LCP (6) is solvable. See the example in Section 5.

######  Remark 2.4.

Assumption 1.1 is also used in [2] for the case $l<m$, which allows more
freedom for the system controls. If $l=m$ and $D$ is invertible, we can write
$Cx(t)+Du(t)-f(t)=-v(t)$ with $v(t)\geq 0$ for a.e. $t\in[0,T]$, where $v\in
L^{2}(0,T)^{l}$. Then the solvability of (5) becomes to find a solution pair
$(x,v)\in H^{1}(0,T)^{n}\times L^{2}(0,T)^{l}$ with $v(t)\geq 0$ satisfying

(8) $\left\\{\begin{aligned}
&\dot{x}(t)=(A-BD^{-1}C)x(t)+BD^{-1}f(t)-BD^{-1}v(t),&\\\
&x(0)=x_{0},\,\,\,x(T)=b,&\\\ \end{aligned}\right.\,\,a.e.\,\,t\in(0,T).$

It then requires the concept of positive controllability [3, 25]. Therefore,
the conditions to ensure the solution set of (8) is nonempty for any
$b\in\mathbb{R}^{n}$ are (i) the matrix

$[BD^{-1}\,\,\,\,(A-BD^{-1}C)BD^{-1}\,\,\,(A-BD^{-1}C)^{2}BD^{-1}\,\,\,\cdot\cdot\cdot\,\,\,(A-BD^{-1}C)^{n-1}BD^{-1}]$

possesses full row rank $n$, and (ii) there is no real eigenvector
$\textbf{w}\in\mathbb{R}^{n}$ of $(A-BD^{-1}C)^{\top}$ such that
$\textbf{w}^{\top}BD^{-1}\textbf{v}\geq 0$ for any
$\textbf{v}\in\mathbb{R}^{m}_{+}$. Then there is a finite time $T_{0}$ such
that the solution set of (5) is nonempty for any $b\in\mathbb{R}^{n}$ and
$T\geq T_{0}$. Hence we can replace Assumption 1.1 in Theorem 2.1 by these two
conditions for the case that $l=m$ and $D$ is invertible.

##  3 Relaxation and sample average approximation (SAA)

In this section, we apply the relaxation and the SAA approach to solve (1). We
consider an independent identically distributed (i.i.d) sample of
$\xi(\omega)$, which is denoted by $\\{\xi_{1},\cdot\cdot\cdot,\xi_{\nu}\\}$,
and use the following relaxation and SAA problem to approximate problem (1):

(9)
$\displaystyle\min_{x,u}\,\,\Phi^{\nu}(x,u):=\frac{1}{\nu}\sum_{\ell=1}^{\nu}F(x(T),\xi_{\ell})+\frac{1}{2}\|x-x_{d}\|_{L^{2}}^{2}+\frac{\delta}{2}\|u-u_{d}\|_{L^{2}}^{2}$
$\displaystyle\textrm{s.t.}\left\\{\begin{aligned} &\left.\begin{aligned}
&\dot{x}(t)=Ax(t)+Bu(t),&\\\ &Cx(t)+Du(t)-f(t)\leq
0,\end{aligned}\right\\}a.e.\,\,t\in(0,T),&\\\
&\left\|\min\left\\{x(T),\frac{1}{\nu}\sum_{\ell=1}^{\nu}[M(\xi_{\ell})x(T)+q(\xi_{\ell})]\right\\}\right\|\leq\epsilon,&\\\
&x(0)=x_{0},\,\,\,\frac{1}{\nu}\sum_{\ell=1}^{\nu}g(x(T),\xi_{\ell})\in
K^{\epsilon}:=\left\\{z\,|\,{\rm
dist}(z,K)\leq\epsilon\right\\},\end{aligned}\right.$

where $\epsilon>0$ is a sufficiently small number.

###  3.1 Convergence of the relaxation and SAA

In this subsection, we show the existence of a solution of problem (9), and
its convergence as $\epsilon\downarrow 0$ and $\nu\to\infty$.

######  Theorem 3.1.

Suppose that the conditions in Theorem 2.1 hold. Then for any $\epsilon>0$,
the SAA problem (9) has an optimal solution
$(x^{\epsilon,\nu},u^{\epsilon,\nu})\in H^{1}(0,T)^{n}\times L^{2}(0,T)^{m}$
for sufficiently large $\nu$ w.p.1.

######  Proof 3.2.

Since the solution set of the linear control system (5) is nonempty for any
$b\in\mathbb{R}^{n}$, for the existence of a feasible solution to the SAA
problem (9), it suffices to show the solution set of the following system

(10) $\left\\{\begin{aligned}
&\left\|\min\left\\{x(T),\frac{1}{\nu}\sum_{\ell=1}^{\nu}[M(\xi_{\ell})x(T)+q(\xi_{\ell})]\right\\}\right\|\leq\epsilon,&\\\
&\frac{1}{\nu}\sum_{\ell=1}^{\nu}g(x(T),\xi_{\ell})\in
K^{\epsilon}\end{aligned}\right.$

is nonempty w.p.1 for any positive number $\epsilon$.

Let $x^{*}(T)$ be a sparse solution of the LCP in (6). From Theorem 2.1, we
know that $x^{*}(T)$ satisfies (6). By the Law of Large Number, for
sufficiently lager $\nu$, $x^{*}(T)$ is a solution of (10). It concludes that
the solution set of the system (10) is nonempty w.p.1 for sufficiently large
$\nu$.

Since $\mathbb{E}[F(\cdot,\xi)]$ is bounded from below, we can also obtain
that $\frac{1}{\nu}\sum_{\ell=1}^{\nu}F(\cdot,\xi_{\ell})$ is bounded from
below with sufficiently large $\nu$. The existence of optimal solutions to
problem (9) is similar to the proof of Theorem 2.1.

Since $F(\cdot,\xi)$ is a smoothing function for a.e. $\xi\in\Xi$, following
the proof of Lemma 3.5 in [7], we can have the following lemma.

######  Lemma 3.3.

Let $\mathcal{C}_{1}\times\mathcal{C}_{2}$ denote a compact subset of
$H^{1}(0,T)^{n}\times L^{2}(0,T)^{m}$. It holds that $\Phi^{\nu}$ epiconverges
to $\Phi$ w.p.1 over $\mathcal{C}_{1}\times\mathcal{C}_{2}$ as
$\nu\rightarrow\infty$.

Let $\mathcal{Z}^{\epsilon,\nu}$ and $\mathcal{Z}$ denote the solution sets of
(10) and (6), respectively. Let $\mathcal{S}^{\epsilon,\nu}$ and $\mathcal{S}$
be the feasible solution sets, and $\hat{\mathcal{S}}^{\epsilon,\nu}$ and
$\hat{\mathcal{S}}$ be optimal solution sets of (9) and (1), respectively.

######  Theorem 3.4.

Suppose that the conditions of Theorem 2.1 hold, then it holds that
$\lim_{\epsilon\downarrow
0}\lim_{\nu\rightarrow\infty}\mathbb{D}(\mathcal{Z}^{\epsilon,\nu},\mathcal{Z})=0$
w.p.1, $\lim_{\epsilon\downarrow
0}\lim_{\nu\rightarrow\infty}\mathbb{D}(\mathcal{S}^{\epsilon,\nu},\mathcal{S})=0$
w.p.1. and $\lim_{\epsilon\downarrow
0}\lim_{\nu\rightarrow\infty}\mathbb{D}(\hat{\mathcal{S}}^{\epsilon,\nu},\hat{\mathcal{S}})=0$
w.p.1.

######  Proof 3.5.

From Theorem 2.1 and Theorem 3.1, $\mathcal{Z}$ and
$\mathcal{Z}^{\epsilon,\nu}$ with any $\epsilon>0$ and large enough $\nu$ are
nonempty. By (7) and
$\mathcal{B}(0,\beta\mathbb{E}[\kappa_{2}(\xi)]\|\mathbb{E}[q(\xi)]\|_{1})\subseteq
K$ in Theorem 2.1, we know that there is a compact set $\mathcal{X}$ such that
$\mathcal{Z}\subseteq\mathcal{X}$ and
$\mathcal{Z}^{\epsilon,\nu}\subseteq\mathcal{X}$ for sufficiently small
$\epsilon$ and sufficiently large $\nu$. Let

$\phi(x(T)):=\min\\{x(T),\mathbb{E}[M(\xi)x(T)+q(\xi)]\\}\quad{\rm
and}\quad\psi(x(T)):=\mathbb{E}[g(x(T),\xi)].$

For $x(T)\in\mathcal{Z}$, $\phi(x(T))=0$ and $\psi(x(T))\in K.$ From (10), for
$x(T)\in\mathcal{Z}^{\epsilon,\nu}$, there are
$v^{\nu}\in\mathbb{R}^{n},w^{\nu}\in\mathbb{R}^{\ell}$ with
$\|v^{\nu}\|\leq\epsilon$ and $\|w^{\nu}\|\leq\epsilon$ such that

$\displaystyle\phi^{\nu}_{\epsilon}(x(T)):=\min\left\\{x(T),\frac{1}{\nu}\sum_{\ell=1}^{\nu}[M(\xi_{\ell})x(T)+q(\xi_{\ell})]\right\\}+v^{\nu}=0,$
$\displaystyle\psi^{\nu}_{\epsilon}(x(T)):=\frac{1}{\nu}\sum_{\ell=1}^{\nu}g(x(T),\xi_{\ell})+w^{\nu}\in
K.$

Since $\phi$ and $\psi$ are continuous, and $M(\cdot),q(\cdot)$ and
$g(x(T),\cdot)$ satisfy Assumption 1.1, we have $\phi^{\nu}_{\epsilon}$ and
$\psi^{\nu}_{\epsilon}$ converge to $\phi$ and $\psi$ uniformly, respectively
on the compact set $\mathcal{X}$, that is,

$\lim_{\epsilon\downarrow
0}\lim_{\nu\to\infty}\max_{x(T)\in\mathcal{X}}\|\phi^{\nu}_{\epsilon}(x(T))-\phi(x(T))\|=0,$

and

$\lim_{\epsilon\downarrow
0}\lim_{\nu\to\infty}\max_{x(T)\in\mathcal{X}}\|\psi^{\nu}_{\epsilon}(x(T))-\psi(x(T))\|=0.$

Therefore, following Theorem 5.12 in [23], $\lim_{\epsilon\downarrow
0}\lim_{\nu\rightarrow\infty}\mathbb{D}(\mathcal{Z}^{\epsilon,\nu},\mathcal{Z})=0$
w.p.1.

Now we show $\lim_{\epsilon\downarrow
0}\lim_{\nu\rightarrow\infty}\mathbb{D}(\mathcal{S}^{\epsilon,\nu},\mathcal{S})=0$
holds w.p.1. Note that $\mathcal{S}^{\epsilon,\nu}$ and $\mathcal{S}$ are two
nonempty closed sets. Obviously, two nonempty closed sets $\mathcal{S}$ and
$\mathcal{S}^{\epsilon,\nu}$ are the solution sets of problem (5) with
terminal sets $\mathcal{Z}$ and $\mathcal{Z}^{\epsilon,\nu}$, respectively.
For any $p\in L^{2}(0,T)^{l}$, the pair $(\|x\|_{H^{1}},\|u\|_{L^{2}})$, where
$(x,u)$ is a solution of problem (4), is uniquely defined by the terminal
point $x(T)$. In addition, it is clear that a solution $(x,u)$ of problem (4)
is continuous with respect to the terminal point $x(T)$. Hence, for any
$(x^{\epsilon,\nu},u^{\epsilon,\nu})\in\mathcal{S}^{\epsilon,\nu}$ and
$(x,u)\in\mathcal{S}$, we have
$(x^{\epsilon,\nu},u^{\epsilon,\nu})\rightarrow(x,u)$ w.p.1 in the norm
$\|\cdot\|_{H^{1}}\times\|\cdot\|_{L^{2}}$ when
$x^{\epsilon,\nu}(T)\rightarrow x(T)$ w.p.1 as $\epsilon\downarrow 0$ and
$\nu\rightarrow\infty$. It then concludes $\lim_{\epsilon\downarrow
0}\lim_{\nu\rightarrow\infty}\mathbb{D}(\mathcal{S}^{\epsilon,\nu},\mathcal{S})=0$
w.p.1.

It is clear that, $\arg\min\Phi=\hat{\mathcal{S}}\subseteq\mathcal{S}$ and
$\arg\min\Phi^{\nu}=\hat{\mathcal{S}}^{\epsilon,\nu}\subseteq\mathcal{S}^{\epsilon,\nu}$
for sufficiently small $\epsilon$ and sufficiently large $\nu$. By
$\lim_{\epsilon\downarrow
0}\lim_{\nu\rightarrow\infty}\mathbb{D}(\mathcal{S}^{\epsilon,\nu},\mathcal{S})=0$
w.p.1, we have, for
$(\hat{x}^{\epsilon,\nu},\hat{u}^{\epsilon,\nu})\in\hat{\mathcal{S}}^{\epsilon,\nu}$,
there is $(\hat{x},\hat{u})\in\mathcal{S}$ such that
$(\hat{x}^{\epsilon,\nu},\hat{u}^{\epsilon,\nu})\rightarrow(\hat{x},\hat{u})$
in the norm $\|\cdot\|_{H^{1}}\times\|\cdot\|_{L^{2}}$ as $\epsilon\downarrow
0$ and $\nu\rightarrow\infty$. In addition, according to Theorem 2.5 in [1],
we obtain $(\hat{x},\hat{u})\in\hat{\mathcal{S}}$ by
$\Phi^{\nu}\xrightarrow{epi}\Phi$ w.p.1, which implies
$\lim_{\epsilon\downarrow
0}\lim_{\nu\rightarrow\infty}\mathbb{D}(\hat{\mathcal{S}}^{\epsilon,\nu},\hat{\mathcal{S}})=0$
w.p.1.

###  3.2 Asymptotics of the SAA optimal value

We introduce the relaxation of problem (1) with a parameter $\epsilon>0$ as
follows

(11) $\displaystyle\min_{x,u}\,\,\Phi(x,u)$
$\displaystyle\textrm{s.t.}\left\\{\begin{aligned} &\left.\begin{aligned}
&\dot{x}(t)=Ax(t)+Bu(t),&\\\ &Cx(t)+Du(t)-f(t)\leq 0,&\\\
\end{aligned}\right\\}a.e.\,\,t\in(0,T),&\\\
&\|\min\\{x(T),\mathbb{E}[M(\xi)x(T)+q(\xi)]\\}\|\leq\epsilon,&\\\
&x(0)=x_{0},\,\mathbb{E}[g(x(T),\xi)]\in K^{\epsilon}.\end{aligned}\right.$

Denote by $\mathcal{Z}^{\epsilon}$ the solution set of the terminal
constraints of (11). Denote by $\mathcal{S}^{\epsilon}$ and
$\hat{\mathcal{S}}^{\epsilon}$ the feasible solution set and optimal solution
set of (11), respectively. Recall that $\mathcal{Z}$ is the solution set of
(6), and $\mathcal{S}$ and $\hat{\mathcal{S}}$ are the feasible solution set
and optimal solution set of (1), respectively. It is clear that
$\mathcal{Z}\subseteq\mathcal{Z}^{\epsilon}$ and
$\mathcal{S}\subseteq\mathcal{S}^{\epsilon}$, which mean that
$\Phi(\hat{x}^{\epsilon},\hat{u}^{\epsilon})\leq\Phi(\hat{x},\hat{u})$ for any
$(\hat{x}^{\epsilon},\hat{u}^{\epsilon})\in\hat{\mathcal{S}}^{\epsilon}$ and
$(\hat{x},\hat{u})\in\hat{\mathcal{S}}$. Therefore, $\mathcal{Z}^{\epsilon}$,
$\mathcal{S}^{\epsilon}$ and $\hat{\mathcal{S}}^{\epsilon}$ are nonempty since
$\mathcal{Z}$ and $\mathcal{S}$ are nonempty.

According to Theorem 2.1 and Remark 2.3, we also conclude that
$\mathcal{Z}^{\epsilon}$ and $\hat{\mathcal{S}}^{\epsilon}$ are compact. It
can also be derived that $\lim_{\epsilon\downarrow
0}\mathbb{D}(\mathcal{Z}^{\epsilon},\mathcal{Z})=0$, $\lim_{\epsilon\downarrow
0}\mathbb{D}(\mathcal{S}^{\epsilon},\mathcal{S})=0$ and
$\lim_{\epsilon\downarrow
0}\mathbb{D}(\hat{\mathcal{S}}^{\epsilon},\hat{\mathcal{S}})=0$. It is clear
that (9) is the corresponding SAA problem of (11). By Theorem 3.4, we conclude
that
$\lim_{\nu\rightarrow\infty}\mathbb{D}({\mathcal{Z}}^{\epsilon,\nu},{\mathcal{Z}}^{\epsilon})=0$,
$\lim_{\nu\rightarrow\infty}\mathbb{D}({\mathcal{S}}^{\epsilon,\nu},{\mathcal{S}}^{\epsilon})=0$
and
$\lim_{\nu\rightarrow\infty}\mathbb{D}(\hat{\mathcal{S}}^{\epsilon,\nu},\hat{\mathcal{S}}^{\epsilon})=0$.

In the rest of this section, we study the asymptotics of optimal value of the
SAA problem (9) for a fixed $\epsilon>0$.

Since $\min\\{x(T),\mathbb{E}[M(\xi)x(T)+q(\xi)]\\}=0$ and
$\mathbb{E}[g(x(T),\xi)]\in K$ for any $x(T)\in\mathcal{Z}$, we have
$\mathcal{Z}\subseteq\textrm{int}\mathcal{Z}^{\epsilon}$, which means that
$\textrm{int}\mathcal{Z}^{\epsilon}\neq\emptyset$. Let

$\hat{\mathcal{Z}}=\\{x(T):\,(x,u)\in\hat{\mathcal{S}}\\}\quad{\rm
and}\quad\hat{\mathcal{Z}}^{\epsilon}=\\{x(T):\,(x,u)\in\hat{\mathcal{S}}^{\epsilon}\\}.$

Obviously, we have
$\hat{\mathcal{Z}}\subseteq\textrm{int}\mathcal{Z}^{\epsilon}$ and
$\lim_{\epsilon\downarrow
0}\mathbb{D}(\hat{\mathcal{Z}}^{\epsilon},\hat{\mathcal{Z}})=0$.

We give the following assumptions. The set $\hat{\mathcal{Z}}$ is a singleton.

  * (i)

There exists a nonnegative measurable function $\kappa_{1}(\xi)$ with
$\mathbb{E}[\kappa_{1}^{2}(\xi)]<\infty$ such that for any $z_{1}$,
$z_{2}\in\mathbb{R}^{n}$ and $\xi\in\Xi$,

$|F(z_{1},\xi)-F(z_{2},\xi)|\leq\kappa_{1}(\xi)\|z_{1}-z_{2}\|,$

and $\mathbb{E}[F^{2}(z,\xi)]<\infty$ for any $z\in\mathbb{R}^{n}$.

  * (ii)

The function $\mathbb{E}[F(\cdot,\xi)]$ is a strongly convex function, that
is, there is a constant $\mu>0$ such that, for any $z_{1}$,
$z_{2}\in\mathbb{R}^{n}$ and $\tau\in(0,1)$,

$\displaystyle\mathbb{E}[F((1-\tau)z_{1}+\tau
z_{2},\xi)]\leq(1-\tau)\mathbb{E}[F(z_{1},\xi)]+\tau\mathbb{E}[F(z_{2},\xi)]-\frac{\mu\tau(1-\tau)}{2}\|z_{1}-z_{2}\|^{2}.$

######  Theorem 3.6.

Suppose that the conditions of Theorem 2.1, Assumption 3.2 and Assumption 3.2
hold. Let $(\hat{x}^{\epsilon},\hat{u}^{\epsilon})$ and
$(\hat{x}^{\epsilon,\nu},\hat{u}^{\epsilon,\nu})$ be optimal solutions of (11)
and (9), respectively. Then for sufficiently small $\epsilon$ and sufficiently
large $\nu$, we have

$\sqrt{\nu}(\Phi^{\nu}(\hat{x}^{\epsilon,\nu},\hat{u}^{\epsilon,\nu})-\Phi(\hat{x}^{\epsilon},\hat{u}^{\epsilon}))\quad\xrightarrow{D}\quad\mathcal{N}(0,\sigma^{2}(\hat{x}^{\epsilon}(T))),$

where “ $\xrightarrow{D}$ ” denotes convergence in distribution and
$\mathcal{N}(0,\sigma^{2}(\hat{x}^{\epsilon}(T)))$ denotes the norm
distribution with mean 0 and variance
$\sigma^{2}(\hat{x}^{\epsilon}(T)):=\mathbb{V}\textrm{ar}[F(\hat{x}^{\epsilon}(T),\xi)]$.

######  Proof 3.7.

Since $\hat{\mathcal{Z}}$ is a singleton,
$\hat{\mathcal{Z}}\subseteq\textrm{int}\mathcal{Z}^{\epsilon}$ and
$\lim_{\epsilon\downarrow
0}\mathbb{D}(\hat{\mathcal{Z}}^{\epsilon},\hat{\mathcal{Z}})=0$, we have
$\hat{\mathcal{Z}}^{\epsilon}\subseteq\textrm{int}\mathcal{Z}^{\epsilon}$ for
sufficiently small $\epsilon$, which means that there is a convex set
$\mathcal{Z}_{\mathcal{X}}$ such that
$\hat{\mathcal{Z}}^{\epsilon}\subseteq\mathcal{Z}_{\mathcal{X}}\subseteq\mathcal{Z}^{\epsilon}$
for sufficiently small $\epsilon$. We can also obtain that
$\hat{\mathcal{Z}}^{\epsilon}$ is a singleton for sufficiently small
$\epsilon$ under Assumption 3.2(ii). We argue it by contradiction. Suppose
$(\hat{x}^{\epsilon},\hat{u}^{\epsilon})$ and
$(\check{x}^{\epsilon},\check{u}^{\epsilon})$ are two optimal solutions of
(11) with $\hat{x}^{\epsilon}(T)\neq\check{x}^{\epsilon}(T)$. Then
$(x^{\epsilon}_{\tau},u^{\epsilon}_{\tau}):=((1-\tau)\hat{x}^{\epsilon}+\tau\check{x}^{\epsilon},(1-\tau)\hat{u}^{\epsilon}+\tau\check{u}^{\epsilon})$
with $\tau\in(0,1)$ is also a feasible solution of (11), since
$x^{\epsilon}_{\tau}(T)\in\mathcal{Z}_{\mathcal{X}}\subseteq\mathcal{Z}^{\epsilon}$.
Moreover,

$\displaystyle\Phi(x^{\epsilon}_{\tau},u^{\epsilon}_{\tau})\leq(1-\tau)\Phi(\hat{x}^{\epsilon},\hat{u}^{\epsilon})+\tau\Phi(\check{x}^{\epsilon},\check{u}^{\epsilon})-\frac{\mu\tau(1-\tau)}{2}\|\hat{x}^{\epsilon}(T)-\check{x}^{\epsilon}(T)\|^{2},$

which means
$\Phi(x^{\epsilon}_{\tau},u^{\epsilon}_{\tau})<\Phi(\hat{x}^{\epsilon},\hat{u}^{\epsilon})$
since
$\Phi(\hat{x}^{\epsilon},\hat{u}^{\epsilon})=\Phi(\check{x}^{\epsilon},\check{u}^{\epsilon})$
and $\hat{x}^{\epsilon}(T)\neq\check{x}^{\epsilon}(T)$. It contradicts the
assumption that $(\hat{x}^{\epsilon},\hat{u}^{\epsilon})$ is an optimal
solution of (11), and then we know that $\hat{\mathcal{Z}}^{\epsilon}$ is a
singleton for sufficiently small $\epsilon$.

In the following argument, $\epsilon>0$ is a fixed number such that
$\hat{\mathcal{Z}}^{\epsilon}$ is singleton and
$\hat{\mathcal{Z}}^{\epsilon}\subseteq\textrm{int}{\mathcal{Z}}^{\epsilon}$.
Denote
$\hat{\mathcal{Z}}^{\epsilon,\nu}=\\{x(T):\,(x,u)\in\hat{\mathcal{S}}^{\epsilon,\nu}\\}$.
We then obtain that
$\lim_{\nu\rightarrow\infty}\mathbb{D}(\hat{\mathcal{Z}}^{\epsilon,\nu},\hat{\mathcal{Z}}^{\epsilon})=0$
and
$\hat{\mathcal{Z}}^{\epsilon}\subseteq\textrm{int}\mathcal{Z}^{\epsilon,\nu}$
for sufficiently large $\nu$ according to
$\lim_{\nu\rightarrow\infty}\mathbb{D}({\mathcal{Z}}^{\epsilon,\nu},{\mathcal{Z}}^{\epsilon})=0$.
Therefore, there is a
$(\hat{x}^{\epsilon,\nu},\hat{u}^{\epsilon,\nu})\in\hat{\mathcal{S}}^{\epsilon,\nu}$
such that $\hat{x}^{\epsilon,\nu}(T)\in\textrm{int}\mathcal{Z}^{\epsilon,\nu}$
with sufficiently large $\nu$, which implies that, for sufficiently large
$\nu$, there is a compact set $\mathcal{X}$ such that
$\hat{\mathcal{Z}}^{\epsilon}\subseteq\mathcal{X}\subseteq\mathcal{Z}^{\epsilon}$
and
$\hat{x}^{\epsilon,\nu}(T)\in\mathcal{X}\subseteq\mathcal{Z}^{\epsilon,\nu}$.

The solution $(x,u)$ of ODE (5) is continuous with respect to the state
terminal value $x(T)$ and the pair $(\|x\|_{H^{1}},\|u\|_{L^{2}})$ is uniquely
defined by $x(T)$. Therefore, there is a compact set $\mathfrak{X}$ such that
$\hat{\mathcal{S}}^{\epsilon}\subseteq\mathfrak{X}\subseteq{\mathcal{S}}^{\epsilon}$
and $\mathfrak{X}\subseteq{\mathcal{S}}^{\epsilon,\nu}$ with
$\hat{\mathcal{S}}^{\epsilon,\nu}\cap\mathfrak{X}\neq\emptyset$. To derive the
error of approximation for optimal value of (9) to that of (11), it suffices
to investigate the error approximation for optimal value of the following
problem

(12) $\displaystyle\min_{(x,u)\in\mathfrak{X}}\Phi(x,u)$

and its SAA problem

(13) $\displaystyle\min_{(x,u)\in\mathfrak{X}}\Phi^{\nu}(x,u),$

where $\Phi$ and $\Phi^{\nu}$ are defined in (11) and (9), respectively.
Clearly, $\mathfrak{X}\subseteq{\mathcal{S}}^{\epsilon}$ with
$\hat{\mathcal{S}}^{\epsilon}\cap\mathfrak{X}\neq\emptyset$ and
$\mathfrak{X}\subseteq{\mathcal{S}}^{\epsilon,\nu}$ with
$\hat{\mathcal{S}}^{\epsilon,\nu}\cap\mathfrak{X}\neq\emptyset$ mean that an
optimal solution of (12) is an optimal solution of (11), and an optimal
solution of (13) is also an optimal solution of (9). Therefore, according to
Theorem 5.7 in [23], we can obtain that, under Assumption 3.2,

$\displaystyle\sqrt{\nu}(\Phi^{\nu}(\hat{x}^{\epsilon,\nu},\hat{u}^{\epsilon,\nu})-\Phi(\hat{x}^{\epsilon},\hat{u}^{\epsilon}))$
$\displaystyle=\sqrt{\nu}(\frac{1}{\nu}\sum_{\ell=1}^{\nu}F(\hat{x}^{\epsilon,\nu}(T),\xi_{\ell})+\frac{1}{2}\|\hat{x}^{\epsilon,\nu}-x_{d}\|_{L^{2}}^{2}+\frac{\delta}{2}\|\hat{u}^{\epsilon,\nu}-u_{d}\|_{L^{2}}^{2}$
$\displaystyle\quad\quad\quad-\mathbb{E}[F(\hat{x}^{\epsilon}(T),\xi)]-\frac{1}{2}\|\hat{x}^{\epsilon}-x_{d}\|_{L^{2}}^{2}-\frac{\delta}{2}\|\hat{u}^{\epsilon}-u_{d}\|_{L^{2}}^{2})$
$\displaystyle\xrightarrow{D}\inf_{(x,u)\in\hat{\mathcal{S}}^{\epsilon}}Y(x,u),$

where $Y(x,u)$ has a normal distribution with mean 0 and variance
$\mathbb{V}\textrm{ar}[F(x(T),\xi)]$ with
$(x,u)\in\hat{\mathcal{S}}^{\epsilon}$. Since
$\hat{\mathcal{Z}}^{\epsilon}=\\{\hat{x}^{\epsilon}(T)\\}$ is a singleton,
$Y(x,u)$ for any $(x,u)\in\hat{\mathcal{S}}^{\epsilon}$ has the same normal
distribution with mean 0 and variance
$\mathbb{V}\textrm{ar}[F(\hat{x}^{\epsilon}(T),\xi)]$. It then concludes our
desired result.

##  4 The time-stepping method

We now adopt the time-stepping method for solving problem (9) with a fixed
sample $\\{\xi_{1},\ldots,\xi_{\nu}\\}$, which uses a finite-difference
formula to approximate the time derivative $\dot{x}$. It begins with the
division of the time interval $[0,T]$ into $N$ subintervals for a fixed step
size $h={T}/{N}=t_{i+1}-t_{i}$ where $i=0,\cdot\cdot\cdot,N-1$. Starting from
$\mathbf{x}_{0}^{\nu}=x_{0}$, we compute two finite sets of vectors
$\\{\mathbf{x}^{\epsilon,\nu}_{1},\mathbf{x}^{\epsilon,\nu}_{2},\cdot\cdot\cdot,\mathbf{x}^{\epsilon,\nu}_{N}\\}\subset\mathbb{R}^{n}$
and
$\\{\mathbf{u}^{\epsilon,\nu}_{1},\mathbf{u}^{\epsilon,\nu}_{2},\cdot\cdot\cdot,\mathbf{u}^{\epsilon,\nu}_{N}\\}\subset\mathbb{R}^{m}$
in the following manner:

(14)
$\displaystyle\min_{\\{\mathbf{x}_{i},\mathbf{u}_{i}\\}_{i=1}^{N}}\,\,\frac{1}{\nu}\sum_{\ell=1}^{\nu}F(\mathbf{x}_{N},\xi_{\ell})+\frac{h}{2}\sum_{i=1}^{N}\left(\|\mathbf{x}_{i}-x_{d,i}\|^{2}+{\delta}\|\mathbf{u}_{i}-u_{d,i}\|^{2}\right)$
$\displaystyle\textrm{s.t.}\left\\{\begin{aligned} &\left.\begin{aligned}
&\mathbf{x}_{i+1}-\mathbf{x}_{i}=hA\mathbf{x}_{i+1}+hB\mathbf{u}_{i+1},&\\\
&C\mathbf{x}_{i+1}+D\mathbf{u}_{i+1}-f_{i+1}\leq
0,\end{aligned}\right\\}i=0,1,\cdot\cdot\cdot,N-1,&\\\
&\left\|\min\left\\{\mathbf{x}_{N},\frac{1}{\nu}\sum_{\ell=1}^{\nu}[M(\xi_{\ell})\mathbf{x}_{N}+q(\xi_{\ell})]\right\\}\right\|\leq\epsilon,&\\\
&\frac{1}{\nu}\sum_{\ell=1}^{\nu}g(\mathbf{x}_{N},\xi_{\ell})\in
K^{\epsilon},\end{aligned}\right.$

where $\epsilon>0$ is a sufficiently small number, $x_{d,i}=x_{d}(t_{i})$,
$u_{d,i}=u_{d}(t_{i})$ and $f_{i}=f(t_{i})$ for $i\in[N]$.

######  Theorem 4.1.

Suppose that the conditions of Theorem 2.1 hold, then for any $\epsilon>0$,
problem (14) has an optimal solution for sufficiently large $\nu$ and
sufficiently small $h$ w.p.1.

######  Proof 4.2.

Theorem 3.1 has shown that the solution set of (10) with any $\epsilon>0$ is
nonempty for sufficiently large $\nu$ w.p.1. About the existence of feasible
solution to problem (14), it suffices to show that the following problem has a
solution for any $b\in\mathbb{R}^{n}$,

(15) $\left\\{\begin{aligned} &\left.\begin{aligned}
&\mathbf{x}_{i+1}=\mathbf{x}_{i}+hA\mathbf{x}_{i+1}+hB\mathbf{u}_{i+1},&\\\
&C\mathbf{x}_{i+1}+D\mathbf{u}_{i+1}-f_{i+1}\leq
0,\end{aligned}\right\\}i=0,1,\cdot\cdot\cdot,N-1,&\\\
&\mathbf{x}_{0}=x_{0},\quad\quad\mathbf{x}_{N}=b.\end{aligned}\right.$

Firstly, denote $A_{h}=I-h(A-BD^{\dagger}C)$. For sufficiently small $h$,
$A_{h}$ is nonsingular. Similar with the proof of Theorem 2.1, from
$\mathbf{x}_{i+1}=\mathbf{x}_{i}+h(A-BD^{\dagger}C)\mathbf{x}_{i+1}+hBD^{\dagger}p_{i+1}$,
the following iteration with $\mathbf{x}_{0}=x_{0}$,

$\mathbf{x}_{i+1}=A_{h}^{-1}(\mathbf{x}_{i}+hBD^{\dagger}p_{i+1}),\quad
i=0,1,\cdot\cdot\cdot,N-1$

generates a solution $\\{\bar{\mathbf{x}}_{i}\\}_{i=1}^{N}$ of the system with
$\mathbf{x}_{0}=x_{0},$

$\mathbf{x}_{i+1}=\mathbf{x}_{i}+hA\mathbf{x}_{i+1}+hB\mathbf{u}_{i+1},\quad
C\mathbf{x}_{i+1}+D\mathbf{u}_{i+1}=p_{i+1},\quad i=0,1,\cdot\cdot\cdot,N-1,$

for any given $p_{i}\in\mathbb{R}^{l}$, $i=1,\dots,N.$

From Assumption 1.1 and the nonsingularity of $A_{h}$, we know that the matrix

$\mathcal{\tilde{R}}_{d}:=[BY\,\,\,\,A_{h}BY\,\,\cdot\cdot\cdot\,\,A_{h}^{n-1}BY]$

has full row rank $n$. Hence the matrix

$\mathcal{R}_{d}:=[hA_{h}^{-1}BY\,\,\,\,h(A_{h}^{-1})^{2}BY\,\,\cdot\cdot\cdot\,\,h(A_{h}^{-1})^{n}BY]$

has full row rank $n$. According to Theorem 3.1.1 in [15], the system with
$\mathbf{x}_{0}=0$,

$\left\\{\begin{aligned}
&\mathbf{x}_{i+1}=A_{h}^{-1}(\mathbf{x}_{i}+hBY{v}_{i+1}),\,\,i=0,1,\cdot\cdot\cdot,N-1,&\\\
&\mathbf{x}_{N}=b-\bar{\mathbf{x}}_{N},\end{aligned}\right.$

admits a solution $\\{\tilde{\mathbf{x}}_{i},\tilde{v}_{i}\\}_{i=1}^{N}$ for
any $b\in\mathbb{R}^{n}$. Therefore,
$\\{\tilde{\mathbf{x}}_{i}+\bar{\mathbf{x}}_{i},\tilde{v}_{i}\\}_{i=1}^{N}$ is
a solution of the following equation

$\left\\{\begin{aligned}
&\mathbf{x}_{i+1}=A_{h}^{-1}(\mathbf{x}_{i}+hBYv_{i+1}+hBD^{\dagger}p_{i+1}),\,\,i=0,1,\cdot\cdot\cdot,N-1,&\\\
&\mathbf{x}_{0}=x_{0},\quad\quad\mathbf{x}_{N}=b.\end{aligned}\right.$

Let
$\tilde{\mathbf{u}}_{i}=Y\tilde{v}_{i}+D^{\dagger}(p_{i}-C(\tilde{\mathbf{x}}_{i}+\bar{\mathbf{x}}_{i}))$.
Then it is easy to verify that
$\\{\tilde{\mathbf{x}}_{i}+\bar{\mathbf{x}}_{i},\tilde{\mathbf{u}}_{i}\\}_{i=1}^{N}$
is a solution of (15) by setting $p_{i}=f_{i}+\tilde{p}_{i}$ for any
$\tilde{p}_{i}\leq 0$.

Since $\mathbb{E}[F(\cdot,\xi)]$ is bounded from below, we can also obtain
that $\frac{1}{\nu}\sum_{\ell=1}^{\nu}F(\cdot,\xi_{\ell})$ is also bounded
from below with sufficiently large $\nu$. Similar to Theorem 5.1 in [2], we
can prove a minimizing sequence tends to an optimal solution of (14), which
shows the existence of optimal solutions to (14) with any $\epsilon>0$ for
sufficiently large $\nu$ and sufficiently small $h$.

Let
$\\{\mathbf{x}^{\epsilon,\nu}_{i},\mathbf{u}^{\epsilon,\nu}_{i}\\}_{i=1}^{N}$
be a solution of (14). We define a piecewise linear function
${x}^{\epsilon,\nu}_{h}$ and a piecewise constant function
${u}^{\epsilon,\nu}_{h}$ on $[0,T]$ as below:

(16)
$\displaystyle{x}^{\epsilon,\nu}_{h}(t)=\mathbf{x}^{\epsilon,\nu}_{i}+\frac{t-t_{i}}{h}(\mathbf{x}^{\epsilon,\nu}_{i+1}-\mathbf{x}^{\epsilon,\nu}_{i}),\quad{u}^{\epsilon,\nu}_{h}(t)=\mathbf{u}^{\epsilon,\nu}_{i+1},\,\,\,\,\forall\,\,t\in(t_{i},t_{i+1}].$

Denote $\hat{\mathcal{S}}^{\epsilon,\nu}_{h}$ the set of
$(\hat{x}_{h}^{\epsilon,\nu},\hat{u}_{h}^{\epsilon,\nu})\in
H^{1}(0,T)^{n}\times L^{2}(0,T)^{m}$, where
$(\hat{x}^{\epsilon,\nu}_{h},\hat{u}^{\epsilon,\nu}_{h})$ are defined in (16)
based on an optimal solution
$\\{\hat{\mathbf{x}}^{\epsilon,\nu}_{i},\hat{\mathbf{u}}^{\epsilon,\nu}_{i}\\}_{i=1}^{N}$
of (14). Let

$\Phi^{\nu}_{h}(x^{\epsilon,\nu}_{h},u^{\epsilon,\nu}_{h})=\frac{1}{\nu}\sum_{\ell=1}^{\nu}F(\mathbf{x}^{\epsilon,\nu}_{N},\xi_{\ell})+\frac{h}{2}\sum_{i=1}^{N}\left(\|\mathbf{x}^{\epsilon,\nu}_{i}-x_{d,i}\|^{2}+{\delta}\|\mathbf{u}^{\epsilon,\nu}_{i}-u_{d,i}\|^{2}\right).$

######  Theorem 4.3.

Suppose that the conditions of Theorem 2.1 hold, then we have

$\lim_{\epsilon\downarrow 0}\lim_{\nu\rightarrow\infty}\lim_{h\downarrow
0}\mathbb{D}(\hat{\mathcal{S}}^{\epsilon,\nu}_{h},\hat{\mathcal{S}})=0,\,\,\textrm{
w.p.1.}$

######  Proof 4.4.

Following [12, 18, 21] and the references therein, we can obtain that
$\Phi^{\nu}_{h}$ epiconverges to $\Phi^{\nu}$ as $h\downarrow 0$ over a
compact subset of $H^{1}(0,T)^{n}\times L^{2}(0,T)^{m}$.

Let
$\\{\hat{\mathbf{x}}^{\epsilon,\nu}_{i},\hat{\mathbf{u}}^{\epsilon,\nu}_{i}\\}_{i=1}^{N}$
be an optimal solution of (14), which means the boundedness of
$\\{\hat{u}^{\epsilon,\nu}_{h_{k}}\\}_{k=1}^{\infty}\subseteq L^{2}(0,T)^{m}$.
Since $L^{2}(0,T)^{m}$ is reflexive, there is a subsequence of
$\\{\hat{u}^{\epsilon,\nu}_{h_{k}}\\}$, which we may assume without loss of
generality to be $\\{\hat{u}^{\epsilon,\nu}_{h_{k}}\\}$ itself, having a weak
limit $\hat{u}^{\epsilon,\nu}_{*}\in L^{2}(0,T)^{m}$. It is easy to see that
$(\hat{x}^{\epsilon,\nu}_{h_{k}},\hat{u}^{\epsilon,\nu}_{h_{k}})$ satisfies
the differential equation
$\dot{x}^{\epsilon,\nu}_{h_{k}}(t)=A\mathbf{x}^{\epsilon,\nu}_{i+1}+Bu^{\epsilon,\nu}_{h_{k}}(t)$
for a.e. $t\in(t_{i},t_{i+1})$ with some $i\in[N]$. Therefore, there is
$\hat{x}^{\epsilon,\nu}_{*}\in H^{1}(0,T)^{n}$ such that
$\hat{x}^{\epsilon,\nu}_{h}\rightarrow\hat{x}^{\epsilon,\nu}_{*}$ in
$H^{1}(0,T)^{n}$ by
$\hat{u}^{\epsilon,\nu}_{h}\rightarrow\hat{u}^{\epsilon,\nu}_{*}$ in
$L^{2}(0,T)^{m}$. According to Theorem 2.5 in [1], we can obtain
$\lim_{h\downarrow
0}\mathbb{D}(\hat{\mathcal{S}}^{\epsilon,\nu}_{h},\hat{\mathcal{S}}^{\epsilon,\nu})=0$
with some $\epsilon>0$ and sufficiently large $\nu$ and then
$\lim_{\epsilon\downarrow 0}\lim_{\nu\rightarrow\infty}\lim_{h\downarrow
0}\mathbb{D}(\hat{\mathcal{S}}^{\epsilon,\nu}_{h},\hat{\mathcal{S}})=0$ w.p.1.

###  4.1 Error estimates of optimal values of problem (14) to problem (9)

In this subsection, we investigate the Euler approximation of problem (9). Our
results are related to the Euler approximation of the optimal control problem
with two-point differential system [13, Theorem 5], which requires the
convexity of the terminal set. However, the terminal constraint set
$\mathcal{Z}^{\epsilon,\nu}$ in (9) is generally nonconvex due to the
existence of the complementarity constraints.

######  Lemma 4.5.

Suppose that the conditions of Theorem 2.1 hold. Let
$({x}^{\epsilon,\nu}_{h},{u}^{\epsilon,\nu}_{h})$ be defined in (16) by a
feasible solution
$\\{\mathbf{x}^{\epsilon,\nu}_{i},\mathbf{u}^{\epsilon,\nu}_{i}\\}_{i=1}^{N}$
of (14). Then, for sufficiently small $h$, there is a feasible solution
$(x^{\epsilon,\nu},u^{\epsilon,\nu})$ of problem (9) such that

$\|x^{\epsilon,\nu}-{x}^{\epsilon,\nu}_{h}\|_{L^{2}}\leq\sqrt{T}C_{y}h,\,\,\,\|x^{\epsilon,\nu}-{x}^{\epsilon,\nu}_{h}\|_{H^{1}}\leq\sqrt{1+T}C_{y}h,\,\,\,\|u^{\epsilon,\nu}-{u}^{\epsilon,\nu}_{h}\|_{L^{2}}\leq
C_{u}h,$

where $C_{y}$ and $C_{u}$ are two constants independent of $h$.

######  Proof 4.6.

We denote two positive constants $\theta_{x}$ and $\theta_{u}$ such that
$\max_{i\in[N]}\|\mathbf{x}^{\epsilon,\nu}_{i}\|\leq\theta_{x}$ and
$\max_{i\in[N]}\|\mathbf{u}^{\epsilon,\nu}_{i}\|\leq\theta_{u}$. According to
Theorem 4.1, there are $v_{i}\in\mathbb{R}^{m-l}$ and $\tilde{p}_{i}\leq 0$
such that
$\mathbf{u}^{\epsilon,\nu}_{i}=Y{v}_{i}+D^{\dagger}(f_{i}+\tilde{p}_{i}-C\mathbf{x}^{\epsilon,\nu}_{i})$
for $i\in[N].$ Let $x^{\epsilon,\nu}(t)$ be the solution of the following
system, for $t\in(t_{i},t_{i+1}]$,

$\left\\{\begin{aligned}
&\dot{x}^{\epsilon,\nu}(t)=(A-BD^{\dagger}C)x^{\epsilon,\nu}(t)+BY({v}_{i+1}+{a_{i+1}(t-t_{i})})+BD^{\dagger}(\tilde{p}_{i+1}+f(t)),&\\\
&{x}^{\epsilon,\nu}(0)=x_{0},\,\,\,x^{\epsilon,\nu}(T)=\mathbf{x}^{\epsilon,\nu}_{N},\end{aligned}\right.$

where $\\{a_{i}\\}_{i=1}^{N}\subset\mathbb{R}^{m-l}$ fulfills

$\displaystyle\mathbf{x}^{\epsilon,\nu}_{N}=e^{(A-BD^{\dagger}C)T}x_{0}+\sum_{i=0}^{N-1}\left[\int_{t_{i}}^{t_{i+1}}e^{(A-BD^{\dagger}C)(T-\tau)}d\tau
B(Yv_{i+1}+D^{\dagger}\tilde{p}_{i+1})\right.$
$\displaystyle\left.+\int_{t_{i}}^{t_{i+1}}e^{(A-BD^{\dagger}C)(T-\tau)}BD^{\dagger}f(\tau)d\tau+\int_{t_{i}}^{t_{i+1}}e^{(A-BD^{\dagger}C)(T-\tau)}BYa_{i+1}(\tau-
t_{i})d\tau\right].$

In addition, we know that $x_{h}^{\epsilon,\nu}$ solves the differential
equation, for any $t\in(t_{i},t_{i+1}]$,

$\displaystyle\left\\{\begin{aligned}
&\dot{x}^{\epsilon,\nu}_{h}(t)=(A-BD^{\dagger}C)x^{\epsilon,\nu}_{h}(t)+BY({v}_{i+1}+a_{i+1}(t-t_{i}))+BD^{\dagger}(\tilde{p}_{i+1}+f(t))+y(t),&\\\
&{x}^{\epsilon,\nu}_{h}(0)=x_{0},\,\,\,x^{\epsilon,\nu}_{h}(T)=\mathbf{x}^{\epsilon,\nu}_{N},\end{aligned}\right.$

where
$y(t)=(A-BD^{\dagger}C)(\mathbf{x}^{\epsilon,\nu}_{i+1}-x_{h}^{\epsilon,\nu}(t))-BYa_{i+1}(t-t_{i})-BD^{\dagger}(f(t)-f_{i+1})$.
Since $f\in L^{2}(0,T)^{l}$, there is a $h_{0}>0$ such that
$\|f(t)-f_{i+1}\|\leq h$ for any $h\in(0,h_{0}]$ and a.e.
$t\in(t_{i},t_{i+1}]$. Let $\tilde{f}(t)=f_{i+1}$ for $t\in(t_{i},t_{i+1}]$,
we then have $\|f-\tilde{f}\|_{L^{2}}\leq\sqrt{T}h$. It means that
$\|y\|_{L^{2}}\leq C_{y}h$ for any $h\in(0,h_{0}]$, where

$C_{y}=\|A-BD^{\dagger}C\|(\|A\|\theta_{x}+\|B\|\theta_{u})\sqrt{T}+\max_{i\in[N]}\|a_{i}\|\|BY\|\sqrt{\frac{T}{3}}+\|BD^{\dagger}\|\sqrt{T}.$

Therefore, we have, for any $t\in(t_{i},t_{i+1}]$,

$\|x^{\epsilon,\nu}-{x}^{\epsilon,\nu}_{h}\|^{2}_{L^{2}}\leq\int_{0}^{T}\int_{0}^{t}\|\dot{x}^{\epsilon,\nu}(\tau)-\dot{x}^{\epsilon,\nu}_{h}(\tau)\|^{2}d\tau
dt=\int_{0}^{T}\int_{0}^{t}\|y(\tau)\|^{2}d\tau dt\leq\|y\|_{L^{2}}^{2}T.$

Hence, according to the definition of $\|\cdot\|_{H^{1}}$, we obtain that

$\|x^{\epsilon,\nu}-{x}^{\epsilon,\nu}_{h}\|_{H^{1}}\leq\sqrt{1+T}\|y\|_{L^{2}}\leq\sqrt{1+T}C_{y}h.$

Let
$u^{\epsilon,\nu}(t)=Y(v_{i+1}+a_{i+1}(t-t_{i}))+D^{\dagger}(\tilde{p}_{i+1}+f(t)-Cx^{\epsilon,\nu}(t))$
for any $t\in(t_{i},t_{i+1}]$. It is clear that
$u^{\epsilon,\nu}_{h}(t)=\mathbf{u}^{\epsilon,\nu}_{i+1}=Yv_{i+1}+D^{\dagger}(f_{i+1}+\tilde{p}_{i+1}-C\mathbf{x}^{\epsilon,\nu}_{i+1})$
for any $t\in(t_{i},t_{i+1}]$. Then we have
$\|u^{\epsilon,\nu}-u^{\epsilon,\nu}_{h}\|_{L^{2}}\leq C_{u}h,$ where

$C_{u}=\max_{i\in[N]}\|a_{i}\|\|Y\|\sqrt{\frac{T}{3}}+\sqrt{T}\|D^{\dagger}\|+\|D^{\dagger}C\|(C_{y}+\|A\|\theta_{x}+\|B\|\theta_{u})\sqrt{T}.$

Clearly, according to the definition of $u^{\nu}(t)$, we can obtain that for
any $t\in(t_{i},t_{i+1}]$,

$Cx^{\epsilon,\nu}(t)+Du^{\epsilon,\nu}(t)-f(t)=\tilde{p}_{i+1}\leq 0,$

which shows that $(x^{\epsilon,\nu},u^{\epsilon,\nu})$ is a feasible solution
of problem (9).

######  Lemma 4.7.

Suppose that the conditions of Theorem 2.1 hold. Let
$(x^{\epsilon,\nu},u^{\epsilon,\nu})$ be a feasible solution of problem (9)
with $\|x^{\epsilon,\nu}(t)\|\leq\theta^{\prime}_{x}$ and
$\|u^{\epsilon,\nu}(t)\|\leq\theta^{\prime}_{u}$ for a.e. $t\in[0,T]$, where
$\theta^{\prime}_{x}$ and $\theta^{\prime}_{u}$ are two positive constants.
Then, for sufficiently small $h$, there is
$(x^{\epsilon,\nu}_{h},{u}^{\epsilon,\nu}_{h})$ defined in (16) by a feasible
solution
$\\{\mathbf{x}^{\epsilon,\nu}_{i},\mathbf{u}^{\epsilon,\nu}_{i}\\}_{i=1}^{N}$
of (14), such that

$\|x^{\epsilon,\nu}-{x}^{\epsilon,\nu}_{h}\|_{L^{2}}\leq\sqrt{T}\tilde{C}_{y}h,\,\,\,\|x^{\epsilon,\nu}-{x}^{\epsilon,\nu}_{h}\|_{H^{1}}\leq\sqrt{1+T}\tilde{C}_{y}h,\,\,\,\|u^{\epsilon,\nu}-{u}^{\epsilon,\nu}_{h}\|_{L^{2}}\leq\tilde{C}_{u}h,$

where $\tilde{C}_{y}$ and $\tilde{C}_{u}$ are two positive constants
independent of $h$.

######  Proof 4.8.

Let $(x^{\epsilon,\nu},u^{\epsilon,\nu})\in H^{1}(0,T)^{n}\times
L^{2}(0,T)^{m}$ be a feasible solution of problem (9), then there are $v\in
L^{2}(0,T)^{m-l}$ and $\tilde{p}\in L^{2}(0,T)^{l}$ with $\tilde{p}(t)\leq 0$
for a.e. $t\in[0,T]$ such that
$u^{\epsilon,\nu}(t)=Yv(t)+D^{\dagger}(\tilde{p}(t)+f(t)-Cx^{\epsilon,\nu}(t))$.
In addition, there are $h_{1}>0$ and a piecewise constant function
$\varphi_{v}(t)=\frac{1}{h}\int_{t_{i}}^{t_{i+1}}v(\tau)d\tau:={\varphi}_{i+1}$
for any $t\in(t_{i},t_{i+1}]$ such that $\|v(t)-\varphi_{v}(t)\|\leq h$ for
a.e. $t\in(t_{i},t_{i+1}]$ with $h\in(0,h_{1}]$. There are also $h_{2}>0$ and
a piecewise constant function
$\varphi_{p}(t)=\frac{1}{h}\int_{t_{i}}^{t_{i+1}}\tilde{p}(\tau)d\tau:=\tilde{\varphi}_{i+1}$
for any $t\in(t_{i},t_{i+1}]$ such that $\|\tilde{p}(t)-\varphi_{p}(t)\|\leq
h$ with $h\in(0,h_{2}]$ and $\varphi_{p}(t)\leq 0$ for a.e.
$t\in(t_{i},t_{i+1}]$.

Recall $A_{h}=I-h(A-BD^{\dagger}C)$. For $i=0,1,\cdot\cdot\cdot,N-1$, let
$\mathbf{x}^{\epsilon,\nu}_{0}=x_{0}$ and

$\displaystyle\mathbf{x}^{\epsilon,\nu}_{i+1}=A_{h}^{-1}(\mathbf{x}^{\epsilon,\nu}_{i}+hBY(\varphi_{i+1}+a_{i+1}h)+hBD^{\dagger}(\tilde{\varphi}_{i+1}+f_{i+1})),$

where $\\{a_{i}\\}_{i=1}^{N}\subset\mathbb{R}^{m-l}$ fulfills

$x^{\epsilon,\nu}(T)=A_{h}^{-N}x_{0}+h\sum_{i=0}^{N-1}A_{h}^{-(i+1)}[BY(\varphi_{i+1}+a_{i+1}h)+hBD^{\dagger}(\tilde{\varphi}_{i+1}+f_{i+1})].$

Let
$\mathbf{u}^{\epsilon,\nu}_{i}=Y(\varphi_{i}+a_{i}h)+D^{\dagger}(\tilde{\varphi}_{i}+f_{i}-C\mathbf{x}^{\epsilon,\nu}_{i})$
for $i\in[N]$. Since $\|x^{\epsilon,\nu}(t)\|\leq\theta^{\prime}_{x}$ and
$\|u^{\epsilon,\nu}(t)\|\leq\theta^{\prime}_{u}$ for a.e. $t\in[0,T]$, there
is a partition to $[0,T]$ such that the sequences
$\\{Y\varphi_{i},D^{\dagger}\tilde{\varphi}_{i}\\}_{i=1}^{N}$ and
$\\{\mathbf{x}^{\epsilon,\nu}_{i},\mathbf{u}^{\epsilon,\nu}_{i}\\}_{i=1}^{N}$
are also bounded for any given $N$. We denote that $\tilde{\theta}_{x}$ and
$\tilde{\theta}_{u}$ are two positive constants such that
$\max_{i\in[N]}\|\mathbf{x}^{\epsilon,\nu}_{i}\|\leq\tilde{\theta}_{x}$ and
$\max_{i\in[N]}\|\mathbf{u}^{\epsilon,\nu}_{i}\|\leq\tilde{\theta}_{u}$. It is
clear that $x^{\epsilon,\nu}_{h}$ satisfies,

$\dot{x}^{\epsilon,\nu}_{h}(t)=(A-BD^{\dagger}C)x^{\epsilon,\nu}_{h}(t)+BYv(t)+BD^{\dagger}(\tilde{p}(t)+f(t))+\tilde{y}(t),\,\,\,t\in(t_{i},t_{i+1}],$

where
$\tilde{y}(t)=(A-BD^{\dagger}C)(\mathbf{x}^{\epsilon,\nu}_{i+1}-x^{\epsilon,\nu}_{h}(t))+BY(\varphi_{i+1}+a_{i+1}h-v(t))+BD^{\dagger}(\tilde{\varphi}_{i+1}+f_{i+1}-\tilde{p}(t)-f(t))$.
It means that $\|\tilde{y}\|_{L^{2}}\leq\tilde{C}_{y}h$ for any
$h\in(0,\min\\{h_{0},h_{1},h_{2}\\}]$, where

$\tilde{C}_{y}=\|A-BD^{\dagger}C\|(\|A\|\tilde{\theta}_{x}+\|B\|\tilde{\theta}_{u})\sqrt{T}+\|BY\|(\sqrt{T}+\max_{i\in[N]}\|a_{i}\|)+2\|BD^{\dagger}\|\sqrt{T}.$

Hence
$\|x^{\epsilon,\nu}-{x}^{\epsilon,\nu}_{h}\|_{L^{2}}\leq\sqrt{T}\|\tilde{y}\|_{L^{2}}\leq\sqrt{T}\tilde{C}_{y}h$
and
$\|x^{\epsilon,\nu}-{x}^{\epsilon,\nu}_{h}\|_{H^{1}}\leq\sqrt{1+T}\tilde{C}_{y}h.$
Moreover, we have
$\|u^{\epsilon,\nu}-u^{\epsilon,\nu}_{h}\|_{L^{2}}\leq\tilde{C}_{u}h,$ where

$\tilde{C}_{u}=\|Y\|(\sqrt{T}+\max_{i\in[N]}\|a_{i}\|)+\|C\|\sqrt{T}(\|A\|\tilde{\theta}_{x}+\|B\|\tilde{\theta}_{u}+\tilde{C}_{y})+2\|D^{\dagger}\|\sqrt{T}.$

Obviously, from the definition of $\mathbf{u}^{\epsilon,\nu}_{i}$, we get
$C\mathbf{x}^{\epsilon,\nu}_{i}+D\mathbf{u}^{\epsilon,\nu}_{i}-f_{i}=\tilde{\varphi}_{i}\leq
0$ $(i\in[N])$, which means that
$\\{\mathbf{x}^{\epsilon,\nu}_{i},\mathbf{u}^{\epsilon,\nu}_{i}\\}_{i=1}^{N}$
is a feasible solution of (14).

######  Lemma 4.9.

Suppose that the conditions of Theorem 2.1 hold. Let
$\\{\mathbf{x}^{\epsilon,\nu}_{i},\mathbf{u}^{\epsilon,\nu}_{i}\\}_{i=1}^{N}$
be a feasible solution of (14) with
$\max_{i\in[N]}\|\mathbf{x}^{\epsilon,\nu}_{i}\|\leq\bar{\theta}_{x}$ and
$\max_{i\in[N]}\|\mathbf{u}^{\epsilon,\nu}_{i}\|\leq\bar{\theta}_{u}$, where
$\bar{\theta}_{x}$ and $\bar{\theta}_{u}$ are two positive constants, and let
$({x}^{\epsilon,\nu}_{h},{u}^{\epsilon,\nu}_{h})$ be defined in (16). Then,
for sufficiently small $h$,

$|\Phi^{\nu}(x^{\epsilon,\nu}_{h},u^{\epsilon,\nu}_{h})-\Phi^{\nu}_{h}(x^{\epsilon,\nu}_{h},u^{\epsilon,\nu}_{h})|\leq
C_{T}h,$

where $C_{T}$ is a positive constant independent of $h$.

######  Proof 4.10.

Since
$\\{\mathbf{x}^{\epsilon,\nu}_{i},\mathbf{u}^{\epsilon,\nu}_{i}\\}_{i=1}^{N}$
is a bounded feasible solution of (14),
$\Phi^{\nu}_{h}(x^{\epsilon,\nu}_{h},u^{\epsilon,\nu}_{h})$ is bounded, which
means that there is a $\theta_{o}>0$ such that
$\max_{i\in[N]}\\{\|\mathbf{x}^{\epsilon,\nu}_{i}-x_{d,i}\|,\|\mathbf{u}_{i}^{\epsilon,\nu}-u_{d,i}\|\\}\leq\theta_{o}$.
Therefore, we have
$\Phi^{\nu}(x^{\epsilon,\nu}_{h},u^{\epsilon,\nu}_{h})-\Phi^{\nu}_{h}(x^{\epsilon,\nu}_{h},u^{\epsilon,\nu}_{h})=W_{1}+W_{2}$,
where

$\displaystyle W_{1}$
$\displaystyle=\frac{1}{2}\sum_{i=0}^{N-1}\int_{t_{i}}^{t_{i+1}}\left(\|x^{\epsilon,\nu}_{h}(t)-x_{d}(t)\|^{2}-\|\mathbf{x}^{\epsilon,\nu}_{i+1}-x_{d}(t_{i+1})\|^{2}\right)dt,$
$\displaystyle W_{2}$
$\displaystyle=\frac{\delta}{2}\sum_{i=0}^{N-1}\int_{t_{i}}^{t_{i+1}}\left(\|u^{\epsilon,\nu}_{h}(t)-u_{d}(t)\|^{2}-\|\mathbf{u}^{\epsilon,\nu}_{i+1}-u_{d}(t_{i+1})\|^{2}\right)dt.$

Note that $x_{d}\in L^{2}(0,T)^{n}$ implies that there is $h_{x}>0$ such that
$\|x_{d}(t_{i+1})-x_{d}(t)\|\leq h$ for a.e. $t\in(t_{i},t_{i+1}]$ with
$h\in(0,h_{x}]$. Then we have

$\displaystyle|W_{1}|\leq\frac{1}{2}\sum_{i=0}^{N-1}\int_{t_{i}}^{t_{i+1}}(\|x^{\epsilon,\nu}_{h}(t)-\mathbf{x}^{\epsilon,\nu}_{i+1}\|$
$\displaystyle+\|x_{d}(t_{i+1})-x_{d}(t)\|)\left(\|x^{\epsilon,\nu}_{h}(t)-\mathbf{x}^{\epsilon,\nu}_{i+1}\|+\|x_{d}(t_{i+1})-x_{d}(t)\|+2\|\mathbf{x}^{\epsilon,\nu}_{i+1}-x_{d}(t_{i+1})\|\right)dt$
$\displaystyle\leq\frac{1}{2}\sum_{i=0}^{N-1}(\|A\|\bar{\theta}_{x}+\|B\|\bar{\theta}_{u}+1)((\|A\|\bar{\theta}_{x}+\|B\|\bar{\theta}_{u}+1)h+2\theta_{o})h^{2}$
$\displaystyle\leq\frac{1}{2}(\|A\|\bar{\theta}_{x}+\|B\|\bar{\theta}_{u}+1)((\|A\|\bar{\theta}_{x}+\|B\|\bar{\theta}_{u}+1)h_{x}+2\theta_{o})Th.$

Moreover, $u_{d}\in L^{2}(0,T)^{m}$ implies that there is $h_{u}>0$ such that
$\|u_{d}(t_{i+1})-u_{d}(t)\|\leq h$ for a.e. $t\in(t_{i},t_{i+1}]$ with
$h\in(0,h_{u}]$. Then

$\displaystyle|W_{2}|$ $\displaystyle\leq$
$\displaystyle\frac{\delta}{2}\sum_{i=0}^{N-1}\int_{t_{i}}^{t_{i+1}}\|u_{d}(t_{i+1})-u_{d}(t)\|\left(\|u_{d}(t_{i+1})-u_{d}(t)\|+2\|\mathbf{u}^{\epsilon,\nu}_{i+1}-u_{d}(t_{i+1})\|\right)dt$
$\displaystyle\leq$ $\displaystyle\frac{\delta}{2}(h_{u}+2\theta_{o})Th.$

It then derives our result for $h\in(0,\min\\{h_{x},h_{u}\\}]$ by denoting

$C_{T}=\frac{1}{2}(\|A\|\bar{\theta}_{x}+\|B\|\bar{\theta}_{u}+1)((\|A\|\bar{\theta}_{x}+\|B\|\bar{\theta}_{u}+1)h_{x}+2\theta_{o})T+\frac{\delta}{2}(h_{u}+2\theta_{o})T.$

######  Theorem 4.11.

Suppose that the conditions of Theorem 2.1 hold. Let
$(\hat{x}^{\epsilon,\nu},\hat{u}^{\epsilon,\nu})$ be an optimal solution of
(9), and let $(\hat{x}^{\epsilon,\nu}_{h},\hat{u}^{\epsilon,\nu}_{h})$ be
defined in (16) associated with an optimal solution
$\\{\hat{\mathbf{x}}^{\epsilon,\nu}_{i},\hat{\mathbf{u}}^{\epsilon,\nu}_{i}\\}_{i=1}^{N}$
of (14). Then, for sufficiently small $h$,

(17)
$|\Phi^{\nu}_{h}(\hat{x}^{\epsilon,\nu}_{h},\hat{u}^{\epsilon,\nu}_{h})-\Phi^{\nu}(\hat{x}^{\epsilon,\nu},\hat{u}^{\epsilon,\nu})|\leq\hat{C}_{T}h,$

where $C_{T}$ is a positive constant independent of $h$.

######  Proof 4.12.

Since
$\\{\hat{\mathbf{x}}^{\epsilon,\nu}_{i},\hat{\mathbf{u}}^{\epsilon,\nu}_{i}\\}_{i=1}^{N}$
is an optimal solution of (14), there is $\psi_{0}$ such that
$\max\\{\|\hat{x}^{\epsilon,\nu}_{h}-x_{d}\|_{L^{2}},\|\hat{u}^{\epsilon,\nu}_{h}-u_{d}\|_{L^{2}}\\}\leq\psi_{0}$,
where $(\hat{x}_{h}^{\epsilon,\nu},\hat{u}_{h}^{\epsilon,\nu})$ is defined in
(16) associated with the sequence
$\\{\hat{\mathbf{x}}^{\epsilon,\nu}_{i},\hat{\mathbf{u}}^{\epsilon,\nu}_{i}\\}_{i=1}^{N}$.
Similarly, $(\hat{x}^{\epsilon,\nu},\hat{u}^{\epsilon,\nu})$ is an optimal
solution of (1), which means that there is $\psi_{1}$ such that
$\max\\{\|\hat{x}^{\epsilon,\nu}-x_{d}\|_{L^{2}},\|\hat{u}^{\epsilon,\nu}-u_{d}\|_{L^{2}}\\}\leq\psi_{1}$.

Following Lemma 4.5-Lemma 4.9, there is $\bar{h}>0$ such that for any
$h\in(0,\bar{h}]$ there is $(x^{\epsilon,\nu},u^{\epsilon,\nu})$, which is a
feasible solution of (9) satisfying
$\|x^{\epsilon,\nu}-\hat{x}^{\epsilon,\nu}_{h}\|_{L^{2}}\leq\sqrt{T}C_{y}h$
and $\|u^{\epsilon,\nu}-\hat{u}^{\epsilon,\nu}_{h}\|_{L^{2}}\leq C_{u}h.$
Moreover, according to Lemma 4.7, for any $h\in(0,\bar{h}]$ there is a
$\\{\mathbf{x}^{\epsilon,\nu}_{i},\mathbf{u}^{\epsilon,\nu}_{i}\\}_{i=1}^{N}$,
which is a feasible solution of (14), such that
$\|\hat{x}^{\epsilon,\nu}-{x}^{\epsilon,\nu}_{h}\|_{L^{2}}\leq\sqrt{T}\tilde{C}_{y}h$
and
$\|\hat{u}^{\epsilon,\nu}-{u}^{\epsilon,\nu}_{h}\|_{L^{2}}\leq\tilde{C}_{u}h,$
where $(x^{\epsilon,\nu}_{h},u^{\epsilon,\nu}_{h})$ is defined in (16) based
on the sequence
$\\{\mathbf{x}^{\epsilon,\nu}_{i},\mathbf{u}^{\epsilon,\nu}_{i}\\}_{i=1}^{N}$.

Then we have
$\Phi^{\nu}_{h}(\hat{x}^{\epsilon,\nu}_{h},\hat{u}^{\epsilon,\nu}_{h})\leq\Phi^{\nu}_{h}({x}^{\epsilon,\nu}_{h},{u}^{\epsilon,\nu}_{h})$,
which means

$\displaystyle\Phi^{\nu}_{h}(\hat{x}^{\epsilon,\nu}_{h},\hat{u}^{\epsilon,\nu}_{h})-\Phi^{\nu}(\hat{x}^{\epsilon,\nu},\hat{u}^{\epsilon,\nu})\leq\Phi^{\nu}_{h}({x}^{\epsilon,\nu}_{h},{u}^{\epsilon,\nu}_{h})-\Phi^{\nu}(\hat{x}^{\epsilon,\nu},\hat{u}^{\epsilon,\nu})$
$\displaystyle\leq|\Phi^{\nu}_{h}({x}^{\epsilon,\nu}_{h},{u}^{\epsilon,\nu}_{h})-\Phi^{\nu}({x}^{\epsilon,\nu}_{h},{u}^{\epsilon,\nu}_{h})|+|\Phi^{\nu}({x}^{\epsilon,\nu}_{h},{u}^{\epsilon,\nu}_{h})-\Phi^{\nu}(\hat{x}^{\epsilon,\nu},\hat{u}^{\epsilon,\nu})|.$

Clearly,

$\displaystyle|\Phi^{\nu}({x}^{\epsilon,\nu}_{h},{u}^{\epsilon,\nu}_{h})-\Phi^{\nu}(\hat{x}^{\epsilon,\nu},\hat{u}^{\epsilon,\nu})|$
$\displaystyle\leq$
$\displaystyle\frac{1}{2}\|{x}^{\epsilon,\nu}_{h}-\hat{x}^{\epsilon,\nu}\|_{L^{2}}(\|{x}^{\epsilon,\nu}_{h}-\hat{x}^{\epsilon,\nu}\|_{L^{2}}+2\|\hat{x}^{\epsilon,\nu}-x_{d}\|_{L^{2}})$
$\displaystyle+\frac{\delta}{2}\|{u}^{\epsilon,\nu}_{h}-\hat{u}^{\epsilon,\nu}\|_{L^{2}}(\|{u}^{\epsilon,\nu}_{h}-\hat{u}^{\epsilon,\nu}\|_{L^{2}}+2\|\hat{u}^{\epsilon,\nu}-u_{d}\|_{L^{2}})$
$\displaystyle\leq$
$\displaystyle\left(\frac{1}{2}\sqrt{T}\tilde{C}_{y}(\sqrt{T}\tilde{C}_{y}\bar{h}+2\psi_{1})+\frac{\delta}{2}\tilde{C}_{u}(\tilde{C}_{u}\bar{h}+2\psi_{1})\right)h.$

Hence, according to Lemma 4.9, we get

$\displaystyle\Phi^{\nu}_{h}(\hat{x}^{\epsilon,\nu}_{h},\hat{u}^{\epsilon,\nu}_{h})-\Phi^{\nu}(\hat{x}^{\epsilon,\nu},\hat{u}^{\epsilon,\nu})$
$\displaystyle\leq$
$\displaystyle\left(\frac{1}{2}\sqrt{T}\tilde{C}_{y}(\sqrt{T}\tilde{C}_{y}\bar{h}+2\psi_{1})+\frac{\delta}{2}\tilde{C}_{u}(\tilde{C}_{u}\bar{h}+2\psi_{1})+C_{T}\right)h$
$\displaystyle=:$ $\displaystyle\beta_{1}h.$

From
$\Phi^{\nu}(\hat{x}^{\epsilon,\nu},\hat{u}^{\epsilon,\nu})\leq\Phi^{\nu}({x}^{\epsilon,\nu},{u}^{\epsilon,\nu})$,
we have

$\displaystyle\Phi^{\nu}(\hat{x}^{\epsilon,\nu},\hat{u}^{\epsilon,\nu})-\Phi_{h}^{\nu}(\hat{x}^{\epsilon,\nu}_{h},\hat{u}^{\epsilon,\nu}_{h})\leq\Phi^{\nu}({x}^{\epsilon,\nu},{u}^{\epsilon,\nu})-\Phi_{h}^{\nu}(\hat{x}^{\epsilon,\nu}_{h},\hat{u}^{\epsilon,\nu}_{h})$
$\displaystyle\leq|\Phi^{\nu}({x}^{\epsilon,\nu},{u}^{\epsilon,\nu})-\Phi^{\nu}(\hat{x}^{\epsilon,\nu}_{h},\hat{u}^{\epsilon,\nu}_{h})|+|\Phi^{\nu}(\hat{x}^{\epsilon,\nu}_{h},\hat{u}^{\epsilon,\nu}_{h})-\Phi_{h}^{\nu}(\hat{x}^{\epsilon,\nu}_{h},\hat{u}^{\epsilon,\nu}_{h})|$

and

$\displaystyle|\Phi^{\nu}({x}^{\epsilon,\nu},{u}^{\epsilon,\nu})-\Phi^{\nu}(\hat{x}^{\epsilon,\nu}_{h},\hat{u}^{\epsilon,\nu}_{h})|$
$\displaystyle\leq$
$\displaystyle\frac{1}{2}\|x^{\epsilon,\nu}-\hat{x}^{\epsilon,\nu}_{h}\|_{L^{2}}(\|x^{\epsilon,\nu}-\hat{x}^{\epsilon,\nu}_{h}\|_{L^{2}}+2\|\hat{x}^{\epsilon,\nu}_{h}-x_{d}\|_{L^{2}})$
$\displaystyle+\frac{\delta}{2}\|u^{\epsilon,\nu}-\hat{u}^{\epsilon,\nu}_{h}\|_{L^{2}}(\|u^{\epsilon,\nu}-\hat{u}^{\epsilon,\nu}_{h}\|_{L^{2}}+2\|\hat{u}^{\epsilon,\nu}_{h}-u_{d}\|_{L^{2}})$
$\displaystyle\leq$
$\displaystyle\left(\frac{1}{2}\sqrt{T}{C}_{y}(\sqrt{T}{C}_{y}\bar{h}+2\psi_{0})+\frac{\delta}{2}{C}_{u}({C}_{u}\bar{h}+2\psi_{0})\right)h.$

It means that

$\displaystyle\Phi^{\nu}(\hat{x}^{\epsilon,\nu},\hat{u}^{\epsilon,\nu})-\Phi^{\nu}_{h}(\hat{x}^{\epsilon,\nu}_{h},\hat{u}^{\epsilon,\nu}_{h})$
$\displaystyle\leq$
$\displaystyle\left(\frac{1}{2}\sqrt{T}{C}_{y}(\sqrt{T}{C}_{y}\bar{h}+2\psi_{0})+\frac{\delta}{2}{C}_{u}({C}_{u}\bar{h}+2\psi_{0})+C_{T}\right)h$
$\displaystyle=:$ $\displaystyle\beta_{2}h.$

Hence (17) holds with $\hat{C}_{T}=\max\\{\beta_{1},\beta_{2}\\}.$

##  5 Numerical experiments

We use the following numerical example to illustrate the theoretical results
obtained in this paper.

(18)
$\displaystyle\min_{x,u}\,\,(\mathbb{E}[\xi_{1}^{2}+\xi_{2}]+1)\|x(T)\|^{2}+\frac{1}{2}\left(\|x\|_{L^{2}}^{2}+\|u\|_{L^{2}}^{2}\right)$
$\displaystyle\textrm{s.t.}\left\\{\begin{aligned} &\left.\begin{aligned}
&\dot{x}_{1}(t)=u_{1}(t),&\\\ &\dot{x}_{2}(t)=x_{2}(t)-u_{2}(t),&\\\
&\dot{x}_{3}(t)=u_{3}(t),&\\\ &\dot{x}_{4}(t)=x_{4}(t)-u_{4}(t),&\\\
&x_{1}(t)+u_{2}(t)\leq 0,&\\\ &{x}_{4}(t)+u_{3}(t)\leq 0,&\\\
\end{aligned}\right\\}a.e.\,\,t\in(0,T),&\\\ &x(0)=(1,1,1,1)^{\top},\quad
0\leq x(T)\perp\mathbb{E}[M(\xi)x(T)+q(\xi)]\geq 0,&\\\
&\left(x_{1}(T)+x_{3}(T),(\mathbb{E}[\xi_{1}]+1)(x_{2}(T)+x_{4}(T))\right)^{\top}\in\mathcal{B}(0,\sqrt{6})\subset\mathbb{R}^{2},\end{aligned}\right.$

where

$q(\xi)=\left(\begin{array}[]{c}3+\xi_{2}\\\ \xi_{1}\\\ 1-\xi_{2}\\\
\xi_{1}+1\\\ \end{array}\right)\quad{\rm and}\quad
M(\xi)=\left(\begin{array}[]{cccc}-2-\xi_{1}&0&-\xi_{2}&-\xi_{1}\\\
0&\xi_{2}&-1&0\\\ 0&-\xi_{1}&\xi_{2}&0\\\ \xi_{2}-1&0&0&\xi_{1}\\\
\end{array}\right).$

We set $T=1$, and $\xi_{1}\sim{\cal{N}}(1,0.01)$ and
$\xi_{2}\sim{\cal{U}}(-1,1)$. It is easy to verify that $\mathbb{E}[M(\xi)]$
is a Z-matrix and the controllability matrix in Assumption 1.1

$\mathcal{R}=\left(\begin{array}[]{cccccccc}0&1&0&0&0&0&0&0\\\
0&0&0&1&0&1&0&1\\\ 0&0&1&0&1&0&1&0\\\ -1&0&-1&0&-1&0&-1&0\end{array}\right),$

is full row rank. We can derive that the solution set of the LCP in (18) is

$\\{(0,0,0,0)^{\top},(1,0,0,0)^{\top},(0,1,1,0)^{\top},(1,1,1,0)^{\top}\\}$

and the solution set of the terminal constraints in (18) is

$\\{(0,0,0,0)^{\top},(1,0,0,0)^{\top},(0,1,1,0)^{\top}\\}.$

By computation using Maple for these three terminal vectors, we obtain that
the values of the objective function in (18) with $(0,0,0,0)^{\top}$,
$(1,0,0,0)^{\top}$ and $(0,1,1,0)^{\top}$ are $35.29712213$, $48.71659453$ and
$25.17501124$, respectively, which means that the optimal solution of (18) has
a unique terminal vector $x(T)=(0,1,1,0)^{\top}$. With the terminal vector
$(0,1,1,0)^{\top}$, we derive an optimal solution

$\displaystyle x_{1}^{*}(t)=$
$\displaystyle(-40.3067\sin(at)+0.3685\cos(at))e^{-ct}+(1.3063\sin(at)+0.6315\cos(at))e^{ct}$
$\displaystyle x_{2}^{*}(t)=$
$\displaystyle(17.379\sin(at)+2.4445\cos(at))e^{-ct}+(3.0042\sin(at)-1.4445\cos(at))e^{ct}$
$\displaystyle x_{3}^{*}(t)=$ $\displaystyle
2.0488e^{-1.618t}+1.8734e^{1.618t}-0.2901e^{0.61805t}-2.6321e^{-0.61805t}$
$\displaystyle x_{4}^{*}(t)=$ $\displaystyle
3.315e^{-1.618t}+0.46938e^{1.618t}-1.1578e^{0.61805t}-1.6266e^{-0.61805t}$
$\displaystyle u_{1}^{*}(t)=$
$\displaystyle(51.113\sin(at)-14.198\cos(at))e^{-ct}+(1.4471\sin(at)+1.2488\cos(at))e^{ct}$
$\displaystyle u_{2}^{*}(t)=$
$\displaystyle(40.3067\sin(at)-0.3685\cos(at))e^{-ct}-(1.3063\sin(at)+0.6315\cos(at))e^{ct}$
$\displaystyle{u}^{*}_{3}(t)=$
$\displaystyle-3.315e^{-1.618t}-0.46938e^{1.618t}+1.1578e^{0.61805t}+1.6266e^{-0.61805t}$
$\displaystyle{u}^{*}_{4}(t)=$ $\displaystyle
8.6789e^{-1.618t}-0.2901e^{1.618t}-0.4423e^{0.61805t}-2.6321e^{-0.61805t}$

of problem (18), where $a=0.34066$ and $c=1.2712$.

It is easy to verify that Assumption 1.1, Assumption 3.4 and Assumption 3.5
hold for the functions
$g(x(T),\xi)=(x_{1}(T)+x_{3}(T),(\xi_{1}+1)(x_{2}(T)+x_{4}(T)))^{\top}$ and
$F(x(T),\xi)=(\xi_{1}^{2}+\xi_{2}+1)\|x(T)\|^{2},$ and random matrix $M(\xi)$
and vector $q(\xi)$. Moreover, conditions of Theorem 2.1 hold, since ${\bf
0}\in{\cal V}$, $\mathbb{E}[M(\xi)]$ is a Z-matrix, and
$K=\mathcal{B}(0,\sqrt{6})\subset\mathbb{R}^{2}.$

We apply the relaxation, the SAA scheme and the time-stepping method to
problem (18). We use Matlab built solver _fmincon_ to solve the discrete
approximation problems of problem (18). Setting $\epsilon=0.00001$, for each
pair $(\nu,h)$ with

$\nu\in\\{500,1000,2000,3000,4000\\},\quad
h\in\\{0.008,0.005,0.004,0.002,0.001\\},$

we generate i.i.d. samples
$\Xi^{\nu,k}=\\{\xi^{k}_{1},\ldots,\xi^{k}_{\nu}\\},k=1,\ldots,10000$. We
solve the discrete problem to find a solution
$(x^{\epsilon,\nu}_{h,k},u^{\epsilon,\nu}_{h,k})$ using each of the samples
$\Xi^{\nu,k}$, $k=1,\ldots,10000$. Then we compute the optimal value of the
discrete problem for each $k$

$\Phi^{\nu,k}_{h}(x^{\epsilon,\nu}_{h,k},u^{\epsilon,\nu}_{h,k})=\frac{1}{\nu}\sum^{\nu}_{i=1}F(x^{\epsilon,\nu}_{h,k}(T),\xi^{k}_{i})-\frac{1}{2}(\|x^{\epsilon,\nu}_{h,k}\|^{2}_{L^{2}}+\|u^{\epsilon,\nu}_{h,k}\|^{2}_{L^{2}}).$

The errors between $\Phi(x^{*},u^{*})=25.17501124$ and the optimal value
$\Phi^{\nu}_{h}(x^{\epsilon,\nu}_{h},u^{\epsilon,\nu}_{h})$ are estimated by

$E_{h}^{\epsilon,\nu}=\frac{1}{10000}\sum^{10000}_{k=1}(\Phi(x^{*},u^{*})-\Phi^{\nu,k}_{h}(x^{\epsilon,\nu}_{h,k},u^{\epsilon,\nu}_{h,k}))^{2}.$

The numerical results are shown in FIG. 1, which verify the convergence
results in Sections 3-4.

Figure 1:  Numerical errors between optimal values of (18) and its discrete
problems

##  6 Conclusions

In this paper, we study the optimal control problem with terminal stochastic
linear complementarity constraints (1), and its relaxation-SAA problem (9) and
the relaxation-SAA-time stepping approximation problem (14). We prove the
existence of feasible solutions and optimal solutions to problem (1) in
Theorem 2.1 under the assumption $E[M(\xi)]$ is a Z-matrix or an adequate
matrix. Under the same assumptions of Theorem 2.1, we prove the existence of
feasible solutions and optimal solutions to (9) and (14). We also show the
convergent properties of these two discrete problems (9) and (14) as the
relaxation parameter $\epsilon\downarrow 0$, the sample size
$\nu\rightarrow\infty$ and mesh size $h\downarrow 0$. Moreover, we provide
asymptotics of the SAA optimal value and the error bound of the time-stepping
method. Problem (1) extends optimal control problem with terminal
deterministic linear complementarity constraints in [2] to stochastic
problems. In [2], Benita and Mehlita derived some stationary points and
constraint qualifications under the assumption that the constrained LCP
(2)-(3) is solvable. Theorem 2.1 gives sufficient conditions for the extension
of solutions of (2)-(3).

## References

  * [1] H. Attouch and R. J.-B. Wets, A convergence theory for saddle functions, Trans. Amer. Math. Soc., 280 (1983), pp. 1–41.
  * [2] F. Benita and P. Mehlitz, Optimal control problems with terminal complemenatrity constraints, SIAM J. Optim., 28 (2018), pp. 3079–3104.
  * [3] R. F. Brammer, Controllability in linear autonomous systems with positive controllers, SIAM J. Control Optim., 10 (1972), pp. 339–353.
  * [4] B. Brogliato, Some perspectives on the analysis and control of complementarity systems, IEEE T. Automat. Contr., 48 (2003), pp. 918–935.
  * [5] M. Brokate and P. Krejčí, Optimal control of ODE systems involving a rate independent variational inequality, Discrete Contin. Dyn. Syst. Ser. B, 18 (2013), pp. 331–348.
  * [6] X. Chen, First order conditions for nonsmooth discretized constrained optimal control problems, SIAM J. Control Optim., 42 (2004), pp. 2004–2015.
  * [7] X. Chen, R. J.-B. Wets, and Y. Zhang, Stochastic variational inequalities: residual minimization smoothing sample average approximations, SIAM J. Optim., 22 (2012), pp. 649–673.
  * [8] X. Chen and S. Xiang, Newton iterations in implicit time-stepping scheme for differential linear complementarity systems, Math. Program., 138 (2013), pp. 579–606.
  * [9] X. Chen and S. Xiang, Sparse solutions of linear complementarity problems, Math. Program., 159 (2016), pp. 539–556.
  * [10] C. Clason, Y. Deng, P. Mehlitz, and U. Prüfert, Optimal control problems with control complementarity constraints: existence results, optimality conditions, and a penalty method, Optim. Methods Softw., 35 (2020), pp. 142–170.
  * [11] R. W. Cottle, J.-S. Pang, and R. E. Stone, The Linear Complementarity Problem, Academic Press, Boston, 1992.
  * [12] A. L. Dontchev, Discrete approximations in optimal control, in Nonsmooth Analysis and Geometric Methods in Deterministic Optimal Control, B. S. Mordukhovich and H. J. Sussmann, eds., New York, NY, 1996, Springer New York, pp. 59–80.
  * [13] A. L. Dontchev and W. W. Hager, Euler approximation of the feasible set, Numer. Funct. Anal. Optim., 15 (1994), pp. 245–261.
  * [14] L. Guo and J. J. Ye, Necessary optimality conditions for optimal control problems with equilibrium constraints, SIAM J. Control Optim., 54 (2016), pp. 2710–2733.
  * [15] C. Heij, A. C. M. Ran, and F. van Schagen, Introduction to Mathematical Systems Theory: Discrete Time Linear Systems, Control and Identification (Second Edition), Birkhäuser Cham, Gewerbestrasse, 2021.
  * [16] P. Lambrianides, Q. Gong, and D. Venturi, A new scalable algorithm for computational optimal control under uncertainty, J. Comput. Phys., 420 (2020), p. 109710.
  * [17] P. Mehlitz and G. Wachsmuth, Bilevel Optimal Control: Existence Results and Stationarity Conditions, Springer International Publishing, Cham, 2020, pp. 451–484.
  * [18] B. S. Mordukhovich and T. Pennanen, Epi-convergent discretization of the generalized bolza problem in dynamic optimization, Optim. Lett., 1 (2007), pp. 379–390.
  * [19] C. Phelps, Q. Gong, J. O. Royset, C. Walton, and I. Kaminer, Consistent approximation of a nonlinear optimal control problem with uncertain parameters, Automatica, 50 (2014), pp. 2987–2997.
  * [20] C. Phelps, J. O. Royset, and Q. Gong, Optimal control of uncertain systems using sample average approximations, SIAM J. Control Optim., 54 (2016), pp. 1–29.
  * [21] O. L. Pironneau and E. Polak, Consistent approximations and approximate functions and gradients in optimal control, SIAM J. Control Optim., 41 (2002), pp. 487–510.
  * [22] A. U. Raghunathan, M. S. Diaz, and L. T. Biegler, An MPEC formulation for dynamic optimization of distillation operations, Comput. Chem. Eng., 28 (2004), pp. 2037–2052.
  * [23] A. Shapiro, D. Dentcheva, and A. Ruszczynski, Lectures on Stochastic Programming: Modeling and Theory, SIAM, Philadelphia, 2009.
  * [24] A. Vieira, B. Brogliato, and C. Prieur, Quadratic optimal control of linear complementarity systems: first order necessary conditions and numerical analysis, IEEE T. Automat. Contr., 65 (2020), pp. 1168–1175.
  * [25] H. Yoshida and T. Tanaka, Positive controllability test for continuous-time linear systems, IEEE T. Automat. Contr., 52 (2007), pp. 1685–1689.