Xiao Li<EMAIL_ADDRESS>
Yutong Li<EMAIL_ADDRESS>
Anouck Girard<EMAIL_ADDRESS>
Ilya Kolmanovsky<EMAIL_ADDRESS>
University of Michigan, Ann Arbor, MI, USA

# System-level Safety Guard: Safe Tracking Control through Uncertain Neural
Network Dynamics Models

###### Abstract

The Neural Network (NN), as a black-box function approximator, has been
considered in many control and robotics applications. However, difficulties in
verifying the overall system safety in the presence of uncertainties hinder
the modular deployment of NN in safety-critical systems. In this paper, we
leverage the NNs as predictive models for trajectory tracking of unknown
dynamical systems. We consider controller design in the presence of both
intrinsic uncertainty and uncertainties from other system modules. In this
setting, we formulate the constrained trajectory tracking problem and show
that it can be solved using Mixed-integer Linear Programming (MILP). The
proposed MILP-based solution enjoys a provable safety guarantee for the
overall system, and the approach is empirically demonstrated in robot
navigation and obstacle avoidance through simulations. The demonstration
videos are available at
https://xiaolisean.github.io/publication/2023-11-01-L4DC2024.

###### keywords:

neural networks, system-level safety, uncertainties, trajectory tracking

##  1 Introduction

Robotic and autonomous driving systems are typically structured as a pipeline
of individual modules, which are designed separately to satisfy corresponding
performance requirements and are yet verified at a system level. With recent
advances in machine learning, Neural Networks have been utilized to learn
either individual modules, e.g., localization ([Li et al.(2021)Li, Du, Zeng,
and Jenkins]), mapping ([Roddick and Cipolla(2020)]), path planning ([Barnes
et al.(2017)Barnes, Maddern, and Posner]), or learn the pipeline in an end-to-
end fashion ([Hu et al.(2023)Hu, Yang, Chen, Li, Sima, Zhu, Chai, Du, Lin,
Wang, et al.]). However, the NNs approximate the desired functionalities as
nonlinear mappings from data introducing intrinsic approximation errors. On a
system level, the performance of an NN module can also be affected by
extrinsic uncertainties from other modules. In safety-critical robotics
applications, the consideration of both intrinsic and extrinsic uncertainties
is crucial to the module design, yet vital to securing safety at the system
level.

Research has been dedicated to developing formal safety guarantees for robotic
applications. The methods based on the Control Barrier Function (CBF) have
been investigated to synthesize barrier functions to enforce safety
constraints ([Wang et al.(2018)Wang, Han, and Egerstedt, Dai et al.(2017)Dai,
Gan, Xia, and Zhan]). The CBFs have also been utilized to avoid collisions in
robotic systems with NN-learned dynamics ([Wei and Liu(2021)]) and with an NN
perception module ([Tong et al.(2023)Tong, Dawson, and Fan]). However, the
notion of uncertainties has not been considered in the literature above. In
this work, we consider trajectory tracking controller design leveraging NNs as
predictive models (see Figure 1). Apart from NN prediction errors, we consider
the presence of uncertainties from other modules, which affect the NN
predictions and, subsequently, impact the controller design and the system’s
safety.

[Hewing et al.(2020)Hewing, Wabersich, Menner, and Zeilinger] provides a
review of techniques used to quantify and incorporate model uncertainties in
control and robotic implementations. Probabilistic ([Berkenkamp et
al.(2016)Berkenkamp, Moriconi, Schoellig, and Krause]) and set-membership-
based methods ([Aswani et al.(2013)Aswani, Gonzalez, Sastry, and Tomlin]) have
been considered to handle the uncertainties in nonlinear dynamic models.
Methods have also been investigated to train NN modules with prediction
confidence levels ([Gal and Ghahramani(2016)]). However, the literature has
been focused on the module-level requirements. In this work, we design the
tracking controller, accounting for system-level safety under both intrinsic
and extrinsic uncertainties.

In a similar setting, robust controllers have been developed, that account for
derived error bounds on perception modules ([Dean et al.(2020)Dean, Matni,
Recht, and Ye, Dean et al.(2021)Dean, Taylor, Cosner, Recht, and Ames, Dean
and Recht(2021)]). We assume that uncertainties are bounded within known sets.
Differently from the existing literature, the uncertain inputs of the NN
predictive model yield the activation status of NN neurons uncertain. We focus
on exploiting the structure of NNs in accounting for uncertainty propagation.
Semi-Definite Programming (SDP) ([Hu et al.(2020)Hu, Fazlyab, Morari, and
Pappas]) has been used to assess the robustness of NN-based solutions. Unlike
[Hu et al.(2020)Hu, Fazlyab, Morari, and Pappas], we consider the setting in
which the NN inputs are affected by uncertainties that belong to a decision-
variable-dependent set.

Figure 1: Trajectory tracking through an NN-learned predictive model: The
controller leverages an NN-learned model (with prediction errors) for
predicting the unknown dynamics and tracks the reference trajectories. The NN
predictions depend on the states $x_{k}$ before the sensing module and the
control $u_{k}$ out of the actuator module, which are not directly accessible
by the controller in this pipeline, and are uncertain quantities due to
sensing noises and actuator disturbances.

The contributions of this paper are as follows: 1) We propose an approach
leveraging an NN-learned dynamic model for robust tracking controller design
subject to system-level safety constraints. The safety and dynamics
constraints are informed by uncertainty set propagation through NN and are
handled using MILP. 2) We consider both intrinsic uncertainties from NN
prediction errors and extrinsic uncertainties present in other modules and
establish theoretical safety guarantees at the system level for the proposed
method. 3) The theoretical properties are further demonstrated in simulations
of collision-avoiding navigation of both omnidirectional robots and vehicles.

This paper is organized as follows: In Section 2, we introduce the assumptions
on the actual dynamics of the system as well as the NN learned dynamics, and
we formulate a robust tracking problem. In Section 3, we present our method to
solve the robust tracking problem, using MILP and its theoretical properties.
In Section 4, we use the proposed method, in combination with a Reachability-
Guided RRT algorithm, to navigate a mass-point omnidirectional robot through a
maze filled with obstacles. In Section 5, we leverage a set-theoretical
localization algorithm that provides vehicle state measurements with
uncertainty bounds, and we use our method to navigate a vehicle while avoiding
collisions. Finally, conclusions are given in Section 6 and proofs of the
theoretical properties are presented in Appendix A.

##  2 Problem Formulation

We consider a discrete-time dynamical system represented by
$x_{k+1}=f(x_{k},u_{k})$, where $x_{k}\in\mathcal{X}$ is the state vector and
$\mathcal{X}\subset\mathbb{R}^{n_{x}}$ is the state feasible set;
$u_{k}\in\mathcal{U}$ is the control and
$\mathcal{U}\subset\mathbb{R}^{n_{u}}$ is the control admissible set; and
$f:\mathcal{X}\times\mathcal{U}\rightarrow\mathcal{X}$ is an unknown nonlinear
mapping. In the sequel, we use lowercase letters, e.g., $x,y,u,w$, to
represent vectors, capital letters, e.g., $W$, to define matrices, and
scripted letters, e.g., $\mathcal{X},\mathcal{U}$, to denote sets. We first
discuss the neural network learned dynamic model for controller design that is
subject to both intrinsic and extrinsic uncertainties in Section 2.1. Then, we
introduce the robust tracking problem ensuring safety under uncertainties in
Section 2.2.

###  2.1 Model Preliminaries

As shown in Figure 1, we approximate the dynamics $f$ with a pre-trained
$\ell-$layer fully connected neural network (NN)
$\tilde{f}:\mathcal{X}\times\mathcal{U}\rightarrow\mathcal{X}$ that admits the
following form

$\begin{array}[]{c}z_{i}=\sigma^{(i)}(\hat{z}_{i}),\quad\hat{z}_{i}=W^{(i)}z_{i-1}+b^{(i)},\quad
i=1,\dots,\ell-1,\\\
z_{0}=[x_{k}^{T},u_{k}^{T}]^{T},\quad\tilde{x}_{k+1}=W^{(\ell)}z_{\ell-1}+b^{(\ell)},\end{array}\vspace{-0.5em}$
(1)

where $W^{(i)}\in\mathbb{R}^{n_{i}\times n_{i-1}}$,
$b^{(i)}\in\mathbb{R}^{n_{i}}$ and
$\sigma^{(i)}:\mathbb{R}^{n_{i}}\rightarrow\mathbb{R}^{n_{i}}$ are the weight
matrix, the bias vector and an element-wise nonlinear activation function in
the $i^{th}$ layer, respectively. The NN uses the inputs $x_{k},u_{k}$ to
compute a prediction $\tilde{x}_{k+1}=\tilde{f}(x_{k},u_{k})$ of the actual
state $x_{k+1}$. We consider the overall system to be subject to uncertainties
(see Figure 1) of three types: the NN prediction $\tilde{x}_{k+1}$ of the
actual state $x_{k+1}$ is subject to an unknown prediction error $w^{x}_{k}$,
i.e.,

$\tilde{x}_{k+1}=x_{k+1}+w^{x}_{k},\;w^{x}_{k}\in\mathcal{W}_{x},\vspace{-0.5em}$
(2)

where $\mathcal{W}_{x}\subset\mathbb{R}^{n_{x}}$ is a bounded set; moreover,
the actual states $x_{k}$ and the actual control input $u_{k}$ are both latent
to the controller due to the presence of noises in other modules. We assume
that a state measurement $y_{k}$ of the actual state $x_{k}$ is available from
the sensing module and admits the form

$y_{k}=x_{k}+w^{y}_{k},\;w^{y}_{k}\in\mathcal{W}_{y},\vspace{-0.5em}$ (3)

where $w^{y}_{k}$ is an unknown measurement noise, and
$\mathcal{W}_{y}\subset\mathbb{R}^{n_{x}}$ is a bounded set. Meanwhile, we
design a feedback controller $\pi:\mathcal{X}\times\mathcal{X}\to\mathcal{U}$
to track the reference state $x^{r}$. The actual control signal $u_{k}$ is
subject to an unknown additive actuator disturbance $w^{u}_{k}$, i.e.,

$u_{k}=\tilde{u}_{k}+w^{u}_{k},\;\tilde{u}_{k}=\pi(y_{k},x^{r}),\;w^{u}_{k}\in\mathcal{W}_{u},\vspace{-0.5em}$
(4)

where $\mathcal{W}_{u}\subset\mathbb{R}^{n_{u}}$ is a bounded set.

###  2.2 Robust Constrained Tracking Problem

Assume that there exists an unsafe subset $\mathcal{X}_{u}\subset\mathcal{X}$,
e.g., representing obstacles, that the system should avoid. Given a reference
$x^{r}\in\mathcal{X}_{s}$ in the safe complement subset
$\mathcal{X}_{s}=\mathcal{X}\backslash\mathcal{X}_{u}$ and a state measurement
$y_{k}$ that obeys the assumption in (3), our task is to design a controller
$\pi(y_{k},x^{r})$ to track the reference $x^{r}$ while keeping the system in
the safe subset $\mathcal{X}_{s}$, i.e., ensuring $x_{k+1}\in\mathcal{X}_{s}$
given $x_{k}\in\mathcal{X}_{s}$. These objectives can be accounted for by a
constrained optimization problem,

$\pi(y_{k},x^{r})=\operatorname*{argmin}\limits_{\tilde{u}_{k}}{\Big{(}\max\limits_{x\in\mathcal{F}\big{(}\mathcal{X}_{k},\;\mathcal{U}_{k}(\tilde{u}_{k})\big{)}}\left\lVert
x-x^{r}\right\rVert_{p}\Big{)}}$ (5a) $\text{subject to:
\quad}\tilde{u}_{k}\in\mathcal{U},\;\mathcal{X}_{k}=(y_{k}\ominus\mathcal{W}_{y})\cap\mathcal{X}_{s},\;\mathcal{U}_{k}=\mathcal{U}_{k}(\tilde{u}_{k})=(\tilde{u}_{k}\oplus\mathcal{W}_{u})\cap\mathcal{U},\quad\quad\quad\quad$
(5b)
$\begin{array}[]{c}\tilde{\mathcal{F}}(\mathcal{X}_{k},\mathcal{U}_{k})=\left\\{\tilde{x}_{k+1}\in\mathcal{X}:\;\tilde{x}_{k+1}=\tilde{f}(x_{k},u_{k}),\forall
x_{k}\in\mathcal{X}_{k},u_{k}\in\mathcal{U}_{k}\right\\},\end{array}$ (5c)
$\mathcal{F}(\mathcal{X}_{k},\mathcal{U}_{k})=\tilde{\mathcal{F}}(\mathcal{X}_{k},\mathcal{U}_{k})\oplus\mathcal{W}_{x},\;\mathcal{F}(\mathcal{X}_{k},\mathcal{U}_{k})\subset\mathcal{X}_{s},$
(5d)

where $\tilde{u}$ is the decision variable, $\left\lVert\cdot\right\rVert_{p}$
represents the $\ell_{p}$ norm, $\ominus$ represents the Pontryagin
difference, and $\oplus$ denotes the Minkowski sum. The term
$y_{k}\ominus\mathcal{W}_{y}$ stands for $\\{y_{k}\\}\ominus\mathcal{W}_{y}$
where we omit the curly brackets for simplicity. Given the measurement $y_{k}$
that satisfies (3), the set $\mathcal{X}_{k}$ in constraint (5b) represents an
uncertainty set that contains the actual state $x_{k}$. Similarly, the set
$\mathcal{U}_{k}$ in constraint (5b) depends on the decision variable
$\tilde{u}_{k}$ and contains the actual control $u_{k}$ based on (4). The set-
valued function
$\tilde{\mathcal{F}}:\mathcal{X}\times\mathcal{U}\rightarrow\mathcal{X}$ in
(5c) calculates the one step ahead reachable set of the learned NN dynamic
model. The reachable set $\mathcal{F}(\mathcal{X}_{k},\mathcal{U}_{k})$ in
(5d) is derived from $\tilde{\mathcal{F}}(\mathcal{X}_{k},\mathcal{U}_{k})$
and contains the actual state $x_{k+1}$ based on (2). Meanwhile, the
constraints in (5d) also guarantee the safety of the next state $x_{k+1}$ for
all the possible uncertainties modeled by (2), (3), (4). The objective defined
by (5a) minimizes the maximum distance between the states in the reachable set
$\mathcal{F}(\mathcal{X}_{k},\mathcal{U}_{k})$ and the reference, thereby
steering the system closer to $x^{r}$. Unlike [Hu et al.(2020)Hu, Fazlyab,
Morari, and Pappas], the input set
$\mathcal{U}_{k}=\mathcal{U}_{k}(\tilde{u}_{k})$ is conditioned on the
decision variable $\tilde{u}_{k}$.

###### Remark

We introduce our methodology in the simplest possible setting, i.e., horizon-
one Model Predictive Control with no control cost ([Wei and Liu(2021)]). The
methodology can be extended to consider multiple prediction horizons, by
augmenting (5c) and (5d) with the following constraints:

$\mathcal{F}^{(i+1)}=\tilde{\mathcal{F}}^{(i+1)}\oplus\mathcal{W}_{x},\;\mathcal{F}^{(i+1)}\subset\mathcal{X}_{s},\;\tilde{\mathcal{F}}^{(i+1)}=\tilde{\mathcal{F}}(\mathcal{F}^{(i)},\mathcal{U}_{k}),\;i=1,2,\dots,N-1,\vspace{-0.5em}$

where $N$ is the number of prediction horizons, and
$\mathcal{F}^{(1)}=\mathcal{F}(\mathcal{X}_{k},\mathcal{U}_{k})$. We consider
$p=1$ in the sequel due to its convenient transcription to linear constraints.
We also note that the objective function can use arbitrary $\ell_{p}$ norm and
consider additive control cost.

##  3 Mixed-integer Linear Programming

To simplify the exposition of the approach, we assume the uncertainty sets in
(2), (3), (4) are hyper-cubes defined according to

$\mathcal{W}_{x}=\left\\{w^{x}\in\mathbb{R}^{n_{x}}:\left|w^{x}\right|\leq\epsilon^{x},\;\epsilon^{x}\in\mathbb{R}^{n_{x}},\;\epsilon^{x}\geq
0\right\\},$ (6a)
$\mathcal{W}_{y}=\left\\{w^{y}\in\mathbb{R}^{n_{x}}:\left|w^{y}\right|\leq\epsilon^{y},\;\epsilon^{y}\in\mathbb{R}^{n_{x}},\;\epsilon^{y}\geq
0\right\\},$ (6b)
$\mathcal{W}_{u}=\left\\{w^{u}\in\mathbb{R}^{n_{u}}:\left|w^{u}\right|\leq\epsilon^{u},\;\epsilon^{u}\in\mathbb{R}^{n_{u}},\;\epsilon^{u}\geq
0\right\\}.$ (6c)

Note that, in principle, one can always find a hypercube overbounding a
bounded set. Moreover, we assume the nonlinear functions in the NN model are
ReLU activation functions

$\sigma^{(i)}(x)={\tt
ReLU}(x)=\max\\{0,x\\},\;i=1,\cdots,\ell,\vspace{-0.5em}$ (7)

that are commonly adopted in different NN architectures and have good
empirical performance ([LeCun et al.(2015)LeCun, Bengio, and Hinton]). With a
small abuse of notation, we use $a\leq b$ to denote the element-wise order
between two vectors $a,b\in\mathbb{R}^{n}$. In the sequel, given a lower bound
$a\in\mathbb{R}^{n}$, an upper bound $b\in\mathbb{R}^{n}$ and $a\leq b$, we
use $[a,b]$ to denote a hypercube in $\mathbb{R}^{n}$. We use $I_{n\times n}$
to represent the identity matrix of size $n$. We use $\mymathbb{1}_{n\times
m}$, and $\mymathbb{0}_{n\times m}$ to denote a matrix with all zeros and
ones, respectively, of size $n\times m$. We neglect the subscript $n\times m$
in calculations assuming that the dimensions are appropriate. Meanwhile, we
consider the state feasible set, control admissible set, and unsafe set
according to

$\mathcal{X}=[\underline{x},\overline{x}],\;\underline{x},\overline{x}\in\mathbb{R}^{n_{x}},\;\underline{x}\leq\overline{x},$
(8a)
$\mathcal{U}=[\underline{u},\overline{u}],\;\underline{u},\overline{u}\in\mathbb{R}^{n_{u}},\;\underline{u}\leq\overline{u},$
(8b)
$\mathcal{X}_{u}=\cup_{i}\mathcal{X}_{u}^{(i)},\;\mathcal{X}_{u}^{(i)}=[\underline{x}_{u}^{(i)},\overline{x}_{u}^{(i)}],\;\underline{x}_{u}^{(i)},\overline{x}_{u}^{(i)}\in\mathbb{R}^{n_{x}},\;\underline{x}\leq\underline{x}_{u}^{(i)}\leq\overline{x}_{u}^{(i)}\leq\overline{x}.$
(8c)

The definition of the unsafe set (8c) enables us to tightly over-bound
obstacles of irregular shapes using unions of hypercubes in the optimization
problem. Then, the equations (5b), (5c), (5d) encode the constraints
associated with the input feasibility, NN structural non-linearity, and system
safety, respectively, and are realized using integer decision variables as
described in the respective Sections 3.1, 3.2, and 3.3. We also focus on the
case when $p=1$ in (5a), i.e., the cost function is based on the $\ell_{1}$
norm, which leads to a MILP in Section 3.4. We prove that the formulated
optimization constraints in this section can guarantee safety in Appendix A.

###  3.1 Constraints Embedding Input Feasibility

Given a measurement $y_{k}$ and the decision variable $\tilde{u}_{k}$, the
constraints (5b) define the sets $\mathcal{X}_{k}$, $\mathcal{U}_{k}$ that
guarantee to contain the actual state $x_{k}$ and the control $u_{k}$ while
the decision variable shall be admissible, i.e., $x_{k}\in\mathcal{X}_{k}$,
$u_{k}\in\mathcal{U}_{k}$, $\tilde{u}_{k}\in\mathcal{U}_{k}$. We embed these
feasibility conditions of the NN inputs using linear inequality and equality
constraints in the following proposition.

######  Proposition 3.1.

Given state measurement $y_{k}$ and a decision variable $\tilde{u}_{k}$,
assume that the unknown actual quantities $x_{k}$ and $u_{k}$ obey (3), (4)
with assumptions in (6b), (6c), (8a), (8b). Let the decision variables
$\tilde{u}_{k}\in\mathbb{R}^{n_{u}}$,
$a_{0},b_{0}\in\mathbb{R}^{n_{x}+n_{u}}$, and
$\delta^{a},\delta^{b}\in\mathbb{R}^{n_{u}}$ satisfy the following constraints

$a_{0,1:n_{x}}=\max\\{\underline{x},y_{k}-\epsilon^{y}\\},\quad
b_{0,1:n_{x}}=\min\\{\overline{x},y_{k}+\epsilon^{y}\\},\quad a_{0}\leq
b_{0},$ (9)
$\begin{array}[]{c}\underline{u}\leq\tilde{u}_{k}\leq\overline{u},\quad\delta^{a}_{j},\delta^{b}_{j}\in\\{0,1\\},\;j=1,\dots,n_{u},\vspace{0.3em}\\\
\left\\{\begin{array}[]{l}a_{0,(n_{x}+1):(n_{x}+n_{u})}\geq\underline{u}\\\
a_{0,(n_{x}+1):(n_{x}+n_{u})}\geq\tilde{u}_{k}-\epsilon^{u}\\\
a_{0,(n_{x}+1):(n_{x}+n_{u})}\leq\underline{u}+M(\mymathbb{1}-\delta^{a})\\\
a_{0,(n_{x}+1):(n_{x}+n_{u})}\leq\tilde{u}_{k}-\epsilon^{u}+M\delta^{a}\end{array}\right.,\left\\{\begin{array}[]{l}b_{0,(n_{x}+1):(n_{x}+n_{u})}\leq\overline{u}\\\
b_{0,(n_{x}+1):(n_{x}+n_{u})}\leq\tilde{u}_{k}+\epsilon^{u}\\\
b_{0,(n_{x}+1):(n_{x}+n_{u})}\geq\overline{u}-M(\mymathbb{1}-\delta^{b})\\\
b_{0,(n_{x}+1):(n_{x}+n_{u})}\geq\tilde{u}_{k}+\epsilon^{u}-M\delta^{b}\end{array}\right.,\end{array}$
(10)

where $\delta^{a}_{j}$ denotes the $j$th element in column vector
$\delta^{a}$, $a_{0,m:n}$ represents the vector containing elements between
row $n$ and row $m$ in $a_{0}$, the constant matrix
$M=\operatorname*{Diag}\left(\max\\{\epsilon^{u},\overline{u}-\underline{u}-\epsilon^{u}\\}\right)$,
and $\operatorname*{Diag}(x)\in\mathbb{R}^{n\times n}$ yields a square matrix
with elements of $x$ on the diagonal and zero anywhere else. Then, it is
guaranteed that
$z_{0}=[x_{k}^{T}\;u_{k}^{T}]^{T}\in\left[a_{0},\;b_{0}\right]$.

The constraints in Proposition 3.1 enforce input feasibility, i.e.,
$\tilde{u}_{k}\in\mathcal{U}_{k}$, and provide bounds $a_{0},b_{0}$ to the NN
input subject to extrinsic uncertainties, i.e.,
$a_{0}\leq[x^{T},\;u^{T}]^{T}\leq b_{0}$ for all
$x\in\mathcal{X}_{k},\;u\in\mathcal{U}_{k}$. Specifically, based on assumption
(6b) and (8a), the constraints (9) imply
$\mathcal{X}_{k}\subseteq[a_{0,1:n_{x}},b_{0,1:n_{x}}]$. Based on assumption
(6c) and (8b), the constraints (10) are equivalent to the following
inequalities

$\max\\{\underline{u},\tilde{u}_{k}-\epsilon^{u}\\}=a_{0,(n_{x}+1)\cdots(n_{x}+n_{u})}\leq
b_{0,(n_{x}+1)\cdots(n_{x}+n_{u})}=\min\\{\overline{u},\tilde{u}_{k}+\epsilon^{u}\\},\vspace{-0.5em}$

thereby
$\mathcal{U}_{k}\subseteq[a_{0,(n_{x}+1):(n_{x}+n_{u})},b_{0,(n_{x}+1):(n_{x}+n_{u})}]$.
We use integer variables $\delta^{a},\delta^{b}$ to move the decision variable
$\tilde{u}_{k}$ out of the nonlinear min/max function: $\delta^{a}_{j}=1$
implies the $j$th element of
$\max\\{\underline{u},\tilde{u}_{k}-\epsilon^{u}\\}$ attains the value of the
$j$th element of $\underline{u}$, otherwise $\delta^{a}_{j}=0$;
$\delta^{b}_{j}=1$ indicates the $j$th element of
$\min\\{\overline{u},\tilde{u}_{k}+\epsilon^{u}\\}$ attains the value of the
$j$th element of $\overline{u}$, otherwise $\delta^{b}_{j}=0$. The proof is
presented in A.1.

###  3.2 Constraints Encoding NN Structural Non-Linearity

Using the bounded sets $\mathcal{X}$ in (8a) and $\mathcal{U}$ in (8b), we can
numerically derive lower and upper bounds
$\underline{\hat{z}}_{i},\overline{\hat{z}}_{i}\in\mathbb{R}^{n_{i}},i=1,\dots,\ell$
on the neuron values $\hat{z}_{i}$ and the output $\tilde{x}_{k+1}$ using
interval arithmetic, i.e.,
$\hat{z}_{i}\in[\underline{\hat{z}}_{i},\overline{\hat{z}}_{i}]$ and
$\tilde{x}_{k+1}\in[\underline{\hat{z}}_{\ell},\overline{\hat{z}}_{\ell}]$. In
the sequel, these derived bounds are used to tighten the constraints and limit
the search region for the optimization solver. Given the decision variables
$a_{0},b_{0}$ in Proposition 3.1 as bounds on the NN input, we can explicitly
encode the variable uncertainty set propagation through the NN defined in (5c)
using the following results:

######  Proposition 3.2.

Given $z_{0}=[x_{k}^{T}\;u_{k}^{T}]^{T}\in\left[a_{0},\;b_{0}\right]$,
consider a NN defined by (1) and (7), and let the decision variables
$a_{i-1},b_{i-1}\in\mathbb{R}^{n_{i-1}}$,
$\hat{a}_{i},\hat{b}_{i}\in\mathbb{R}^{n_{i}}$,
$\delta^{--}_{i},\delta^{-+}_{i},\delta^{++}_{i}\in\\{0,1\\}^{n_{i}}$,
$i=1,\dots,\ell-1$, and $a_{k+1},b_{k+1}\in\mathbb{R}^{n_{x}}$ satisfy the
following constraints

$\begin{array}[]{c}\hat{a}_{i,j}=w^{(i)}_{j}S_{i}\left((w^{(i)}_{j})^{T}\right)\left[\begin{array}[]{c}a_{i-1}\\\
b_{i-1}\end{array}\right]+b^{(i)}_{j},\quad\hat{b}_{i,j}=w^{(i)}_{j}S_{i}\left((w^{(i)}_{j})^{T}\right)\left[\begin{array}[]{c}b_{i-1}\\\
a_{i-1}\end{array}\right]+b^{(i)}_{j},\vspace{0.3em}\\\
\underline{\hat{z}}_{i}\leq\hat{a}_{i}\leq\hat{b}_{i}\leq\overline{\hat{z}}_{i},\quad\forall
j=1,\cdots,n_{i},\quad\forall i=1,\dots,\ell-1,\end{array}$ (11)
$\begin{array}[]{c}\left\\{\begin{array}[]{l}a_{i}\geq\hat{a}_{i}\\\
a_{i}\leq\hat{a}_{i}-\operatorname*{Diag}(\underline{\hat{z}}_{i})(\delta^{--}_{i}+\delta^{-+}_{i})\\\
a_{i}\leq\operatorname*{Diag}(\overline{\hat{z}}_{i})\delta^{++}_{i}\end{array}\right.,\quad\left\\{\begin{array}[]{l}b_{i}\geq\hat{b}_{i}\\\
b_{i}\leq\hat{b}_{i}-\operatorname*{Diag}(\underline{\hat{z}}_{i})\delta^{--}_{i}\\\
b_{i}\leq\operatorname*{Diag}(\overline{\hat{z}}_{i})(\delta^{-+}_{i}+\delta^{++}_{i})\end{array}\right.,\vspace{0.3em}\\\
0\leq a_{i}\leq
b_{i},\quad\delta^{--}_{i,j},\delta^{-+}_{i,j},\delta^{++}_{i,j}\in\\{0,1\\},\quad\delta^{--}_{i,j}+\delta^{-+}_{i,j}+\delta^{++}_{i,j}=1,\\\
\forall j=1,\cdots,n_{i},\quad\forall i=1,\dots,\ell-1,\end{array}$ (12)
$\begin{array}[]{c}\underline{\hat{z}}_{\ell}\leq a_{k+1}\leq
b_{k+1}\leq\overline{\hat{z}}_{\ell},\quad
a_{k+1,j}=w^{(\ell)}_{j}S_{\ell}\left((w^{(\ell)}_{j})^{T}\right)\left[\begin{array}[]{c}a_{\ell-1}\\\
b_{\ell-1}\end{array}\right]+b^{(\ell)}_{j},\\\
b_{k+1,j}=w^{(\ell)}_{j}S_{\ell}\left((w^{(\ell)}_{j})^{T}\right)\left[\begin{array}[]{c}b_{\ell-1}\\\
a_{\ell-1}\end{array}\right]+b^{(\ell)}_{j},\quad\forall
j=1,\cdots,n_{i},\end{array}$ (13)

where $\hat{a}_{i,j}$ is the $j^{th}$ element of $\hat{a}_{i}$; $w^{(i)}_{j}$
is the $j^{th}$ row of $W^{(i)}$; $b^{(i)}_{j}$ is the $j^{th}$ element of
$b^{(i)}$; $\delta^{--}_{i,j}$ denotes the $j$th element of $\delta^{--}_{i}$;
$a_{k+1,j},b_{k+1,j}$ represent the $j^{th}$ element of $a_{k+1},b_{k+1}$,
respectively; The functions
$S_{i}:\mathbb{R}^{n_{i-1}}\rightarrow\mathbb{R}^{n_{i-1}\times 2n_{i-1}}$ and
$s:\mathbb{R}\rightarrow\mathbb{R}^{2n_{i-1}}$ are defined according to

$\begin{array}[]{c}S_{i}\left(\left[\begin{array}[]{c}\vdots\\\ w_{q}\\\
\vdots\end{array}\right]\right)=\left[\begin{array}[]{c}\vdots\\\
\left(s(w_{q})\right)^{T}\\\ \vdots\end{array}\right],\quad
s(w_{q}):=\left\\{\begin{array}[]{cc}\left[\begin{array}[]{c}e_{q}\\\
\mymathbb{0}_{n_{i-1}\times 1}\end{array}\right]&\text{if }w_{q}\geq 0\\\
\left[\begin{array}[]{c}\mymathbb{0}_{n_{i-1}\times 1}\\\
e_{q}\end{array}\right]&\text{if
}w_{q}<0\end{array}\right.,\end{array}\vspace{-0.5em}$

and $e_{q}\in\mathbb{R}^{n_{i-1}}$ ($q=1,\dots,n_{i-1}$) has 1 in the $q^{th}$
element and zero anywhere else. Then, the variable reachable set
$\tilde{\mathcal{F}}(\mathcal{X}_{k},\mathcal{U}_{k})$ defined in (5c) is a
subset of hypercube $[a_{k+1},b_{k+1}]$, i.e.,
$\tilde{\mathcal{F}}(\mathcal{X}_{k},\mathcal{U}_{k})\subseteq[a_{k+1},b_{k+1}]$.

From the $(i-1)$th to $i$th layer of the NN, the uncertainty set propagation
is realized through variables $a_{i-1},b_{i-1}$ and $\hat{a}_{i},\hat{b}_{i}$
that are the lower and upper bounds of $z_{i-1}$ and $\hat{z}_{i}$,
respectively (i.e., $z_{i-1}\in[a_{i-1},b_{i-1}]$ and
$\hat{z}_{i}\in[\hat{a}_{i},\hat{b}_{i}]$). The constraints (11) and (13)
encode the uncertainty propagation through the fully connected layers,
$\hat{z}_{i}=W^{(i)}z_{i-1}+b^{(i)}$ and
$\tilde{x}_{k+1}=W^{(\ell)}z_{\ell-1}+b^{(\ell)}$, respectively. The
uncertainty set propagation through the nonlinear ReLU function is enforced in
constraints (12). In constraints (11) and (13), the $q$th row of matrix
$S_{i}((w^{(i)}_{j})^{T})$ switches the upper and lower bounds of $z_{i-1}$ if
the $q$th element in $w^{(i)}_{j}$ is negative. Meanwhile, since the values
$z_{i-1}$ in neurons fall into a bounded hypercube $[a_{i-1},b_{i-1}]$, the
activation status of each ReLU activation function is uncertain. Inspired by
[Tjeng et al.(2017)Tjeng, Xiao, and Tedrake], we introduce integer variables
$\delta^{--}_{i},\delta^{-+}_{i},\delta^{++}_{i}$ in constraints (12) to
encode the uncertainty in the ReLU activation status according to

$\hat{a}_{i,j}\leq\hat{b}_{i,j}\leq 0\text{ if
}\delta^{--}_{i,j}=1;\quad\hat{a}_{i,j}\leq 0\leq\hat{b}_{i,j}\text{ if
}\delta^{-+}_{i,j}=1;\quad 0\leq\hat{a}_{i,j}\leq\hat{b}_{i,j}\text{ if
}\delta^{++}_{i,j}=1.\vspace{-0.5em}$

Furthermore, as can be shown, the hypercube defined by lower bound $a_{k+1}$
and upper bound $b_{k+1}$ is an over-estimation of the actual reachable set
$\tilde{\mathcal{F}}$, i.e., $\tilde{\mathcal{F}}\subseteq[a_{k+1},b_{k+1}]$.
The proof is presented in Appendix A.2.

###  3.3 Constraints Enforcing System Safety

Given a hypercube $[a_{k+1},b_{k+1}]$ as an over-estimation of the reachable
set $\tilde{\mathcal{F}}$ in Proposition 3.2, we enforce the safety constrains
(5d) of the system such that
$x_{k+1}\in\mathcal{F}(\mathcal{X}_{k},\mathcal{U}_{k})\subset\mathcal{X}_{s}$.
We consider the presence of a simple unsafe subset
$\mathcal{X}_{u}=[\underline{x}_{u},\overline{x}_{u}]$ in the following result
and discuss the extension to a union of $\mathcal{X}_{u}^{(i)}$ defined in
(8c) at the end of this section.

######  Proposition 3.3.

Given
$\tilde{\mathcal{F}}(\mathcal{X}_{k},\mathcal{U}_{k})\subseteq[a_{k+1},b_{k+1}]$
and state space defined in (8a), assume that the NN predictions are subject to
bounded additive errors defined in (2) and (6a), and the decision variables
$\underline{x}_{k+1},\overline{x}_{k+1}\in\mathbb{R}^{n_{x}}$,
$\delta^{u}_{1},\delta^{u}_{2}\in\mathbb{R}^{n_{x}}$ satisfy the following
constraints:

$\underline{x}\leq\underline{x}_{k+1}\leq\overline{x}_{k+1}\leq\overline{x},\quad\underline{x}_{k+1}=a_{k+1}-\epsilon^{x},\quad\overline{x}_{k+1}=b_{k+1}+\epsilon^{x},$
(14)
$\begin{array}[]{c}\left\\{\begin{array}[]{l}\overline{x}_{k+1}\leq\overline{x}+\operatorname*{Diag}\left(\underline{x}_{u}-\overline{x}\right)\delta^{u}_{1}\\\
\overline{x}_{k+1}\geq\underline{x}_{u}-\operatorname*{Diag}\left(\underline{x}_{u}-\underline{x}\right)\delta^{u}_{1}\\\
\underline{x}_{k+1}\geq\underline{x}+\operatorname*{Diag}\left(\overline{x}_{u}-\underline{x}\right)\delta^{u}_{2}\\\
\underline{x}_{k+1}\leq\overline{x}_{u}-\operatorname*{Diag}\left(\overline{x}_{u}-\overline{x}\right)\delta^{u}_{2}\\\
\end{array}\right.,\quad\begin{array}[]{c}\delta^{u}_{1,j},\delta^{u}_{2,j}\in\\{0,1\\},\quad\forall
j=1,\dots,n_{x},\\\ \delta^{u}_{1,j}+\delta^{u}_{2,j}\leq 1,\quad\forall
j=1,\dots,n_{x},\\\
\sum\limits_{j=1}^{n_{x}}\left(\delta^{u}_{1,j}+\delta^{u}_{2,j}\right)\geq
1,\end{array}\end{array}$ (15)

where $\delta^{u}_{1,j}$, $\delta^{u}_{2,j}$ are the $j$th element of
$\delta^{u}_{1},\delta^{u}_{2}$, respectively. Then,
$x_{k+1}\in\mathcal{F}(\mathcal{X}_{k},\mathcal{U}_{k})$ where
$\mathcal{F}(\mathcal{X}_{k},\mathcal{U}_{k})$ is the reachable set defined in
(5d), and
$\mathcal{F}(\mathcal{X}_{k},\mathcal{U}_{k})\subseteq[\underline{x}_{k+1},\overline{x}_{k+1}]\subset\mathcal{X}_{s}$
with $\mathcal{X}_{u}=[\underline{x}_{u},\overline{x}_{u}]$.

In Proposition 3.3, the constraints (14) are equivalent to
$[\underline{x}_{k+1},\overline{x}_{k+1}]=[a_{k+1},b_{k+1}]\oplus\mathcal{W}_{x}$
and $[\underline{x}_{k+1},\overline{x}_{k+1}]\subset\mathcal{X}$. Considering
the result
$\tilde{\mathcal{F}}(\mathcal{X}_{k},\mathcal{U}_{k})\subseteq[a_{k+1},b_{k+1}]$
from Proposition 3.2, it is obvious that
$\mathcal{F}(\mathcal{X}_{k},\mathcal{U}_{k})\subseteq[\underline{x}_{k+1},\overline{x}_{k+1}]$
based on the definition of $\mathcal{F}(\mathcal{X}_{k},\mathcal{U}_{k})$ in
(5d). Subsequently, the safety constraint of
$\mathcal{F}(\mathcal{X}_{k},\mathcal{U}_{k})\subset\mathcal{X}_{s}$ can be
enforced with
$[\underline{x}_{k+1},\overline{x}_{k+1}]\subset\mathcal{X}_{s}$, which is
equivalent to $[\underline{x}_{k+1},\overline{x}_{k+1}]\subset\mathcal{X}$ and
$[\underline{x}_{k+1},\overline{x}_{k+1}]\cap\mathcal{X}_{u}=\varnothing$. The
constraints (15) enforce
$[\underline{x}_{k+1},\overline{x}_{k+1}]\cap\mathcal{X}_{u}=\varnothing$
using integer variables according to

$\left\\{\begin{array}[]{cc}\underline{x}_{k+1,j}\leq\overline{x}_{k+1,j}\leq\underline{x}_{u,j}&\text{if
}\delta^{u}_{1,j}=1,\;\delta^{u}_{2,j}=0\\\
\overline{x}_{k+1,j}\geq\underline{x}_{k+1,j}\geq\overline{x}_{u,j}&\text{if
}\delta^{u}_{1,j}=0,\;\delta^{u}_{2,j}=1\\\
\underline{x}_{u,j}\leq\overline{x}_{k+1,j}\leq\overline{x},\;\underline{x}\leq\underline{x}_{k+1,j}\leq\overline{x}_{u,j}&\text{if
}\delta^{u}_{1,j}=0,\;\delta^{u}_{2,j}=0\end{array}\right.,\vspace{-0.5em}$

where
$\overline{x}_{k+1,j},\underline{x}_{k+1,j},\overline{x}_{u,j},\underline{x}_{u,j}$
are the $j^{th}$ element of
$\overline{x}_{k+1},\underline{x}_{k+1},\overline{x}_{u},\underline{x}_{u}$,
respectively. The integer constraints imply that there exists at least one
dimension $j$ where
$[\underline{x}_{k+1,j},\overline{x}_{k+1,j}]\cap[\underline{x}_{u,j},\overline{x}_{u,j}]=\varnothing$.
Subsequently, the hypercube $[\underline{x}_{k+1},\overline{x}_{k+1}]$ has
zero overlaps with the unsafe subset $[\underline{x}_{u},\overline{x}_{u}]$.
Notably, in the case of a complex unsafe region $\mathcal{X}_{u}$, we can
derive a union of hypercubes $\cup_{i}\mathcal{X}_{u}^{(i)}$ as in (8c) that
over-bounds $\mathcal{X}_{u}$. Thereafter, to ensure safety, we can formulate
similar constraints (15) with each individual hypercube
$\mathcal{X}_{u}^{(i)}$ in the union. The proof is available in Appendix A.3.

###  3.4 Safe Tracking Control using MILP

The optimization objective in (5a) is designed to minimize the maximum
distance between the points in the reachable set and the reference state
$x^{r}$. Focusing on the case when $p=1$, i.e., the cost is defined using
$\ell_{1}$ norm, we can introduce a vector of slack variables
$\lambda\in\mathbb{R}^{n_{x}}$, and reformulate (5a) as

$\begin{array}[]{c}\operatorname*{argmin}_{\tilde{u}_{k}}\sum_{q=1}^{n_{x}}\lambda_{q},\end{array}\vspace{-0.5em}$
$\text{subject to:\quad}\begin{array}[]{c}\lambda\geq
0,\;-\lambda\leq\underline{x}_{k+1}-x^{r}\leq\lambda,\;-\lambda\leq\overline{x}_{k+1}-x^{r}\leq\lambda,\end{array}$
(16)

where $\lambda_{q}$ designates the $q$th element of vector $\lambda$. This
objective function relies on the fact that the maximum distance, between a
reference $x^{r}$ and points in the hypercube
$[\underline{x}_{k+1},\overline{x}_{k+1}]$, is attained at the points located
at the boundary of the hypercube. Then, the optimization problem (5) for
tracking the reference state $x^{r}$ safely can be rewritten into a MILP
according to

######  Problem 3.4.

Robust Constrained Tracking Control Problem

$\begin{array}[]{c}\operatorname*{argmin}\limits_{\begin{subarray}{c}\tilde{u}_{k},\;\delta^{a},\;\delta^{b},\;a_{0},\;b_{0},\;a_{k+1},\;b_{k+1},\;\underline{x}_{k+1},\;\overline{x}_{k+1},\;\delta^{u}_{1},\;\delta^{u}_{2},\;\lambda,\\\
a_{i},\;b_{i},\;\hat{a}_{i},\;\hat{b}_{i},\;\delta_{i}^{--},\;\delta_{i}^{-+},\;\delta_{i}^{++},\;i=1,\dots,\ell-1,\end{subarray}}\quad\sum\limits_{q=1}^{n_{x}}\lambda_{q},\\\
\text{subject to: \eqref{eq:mip_input_x}, \eqref{eq:mip_input_u},
\eqref{eq:mip_struct_linear}, \eqref{eq:mip_struct_relu},
\eqref{eq:mip_struct_out}, \eqref{eq:mip_actual_reachable},
\eqref{eq:mip_out_unsafe},
\eqref{eq:mip_objective}.}\end{array}\vspace{-0.5em}$ (17)

We note that the number of decision variables in the MILP problem (17) scales
linearly with the number of neurons in the NN. In the subsequent examples, we
use the YALMIP toolbox ([Lofberg(2004)]) for MATLAB to solve the optimization
(17). It also has the following property. The proof is presented in Appendix
A.

######  Proposition 3.5.

Consider a NN-learned dynamical system defined by (1), (7) that takes control
$u_{k}$ and state $x_{k}$ as inputs and yields $\tilde{x}_{k+1}$ as a
prediction of the next state $x_{k+1}$ and satisfies the bounded additive
error assumption in (2), (6a). We assume that $x_{k}\in\mathcal{X}_{s}$ is
unknown but belongs to a safe set
$\mathcal{X}_{s}=\mathcal{X}\backslash\mathcal{X}_{u}$ that is the complement
set of the unsafe region $\mathcal{X}_{u}$ given by (8c) in the state space
$\mathcal{X}$ defined by (8a). Also assume that a measurement $y_{k}$ of
$x_{k}$ is given that satisfies the assumptions in (3), (6b). If there exists
a solution of the Problem (17) such that the corresponding $\tilde{u}_{k}$ is
in the admissible set $\mathcal{U}$ defined by (8b), then for all actuator
disturbances $w_{k}^{u}$ in set $\mathcal{W}_{u}$ defined by (6c), the actual
control $u_{k}$ subject to this additive disturbance $w_{k}^{u}$ according to
(4) renders the actual next system state safe, i.e.,
$x_{k+1}\in\mathcal{X}_{s}$.

##  4 Obstacle Avoidance and Reachability-Guided RRT

Figure 2: Schematic of obstacle avoidance using an omnidirectional robot:
(Left) The reachability-guided RRT algorithm expands the tree from the start
$x_{0}$ to the goal $x_{g}$ over the safe state space $\mathcal{X}_{s}$. Then,
we use the Dijkstra planning algorithm to find a reference path (black lines
with red dots) from $x_{0}$ to $x_{g}$ that has the shortest distance defined
by $\ell_{1}$ norm. (Middle) We use the proposed method in (17) to track the
reference states (red dots). Our method can guarantee that the robot is
collision-free under uncertainties, i.e., the unknown actual states in blue
asterisks are in the safe set $\mathcal{X}_{s}$, even with most of the
reference points located near the obstacles.

We consider an omnidirectional robot as a mass point and use the following
equation to represent its kinematics:
$x_{k+1}=f(x_{k},u_{k})=x_{k}+u_{k}+w_{k}$, where $x_{k}\in\mathbb{R}^{2}$ is
the robot coordinates in the $X-Y$ plane and $u_{k}\in\mathbb{R}^{2}$ contains
the displacements, and the values of the disturbance $w_{k}\sim
U(-\epsilon^{x},\epsilon^{x})$ are sampled uniformly from the interval
$[-\epsilon^{x},\epsilon^{x}]$. We use the sets
$\mathcal{U}=\left\\{[u_{1}\;u_{2}]^{T}\in\mathbb{R}^{2}:-0.25\leq
u_{1},u_{2}\leq-0.25\right\\}$,
$\mathcal{X}=\left\\{[x_{1}\;x_{2}]^{T}\in\mathbb{R}^{2}:-1\leq
x_{1},x_{2}\leq 10\right\\}$ together with obstacles $\mathcal{X}_{u}^{(i)}$
visualized in orange boxes in Figure 2. The NN is manually constructed and
admits the following form

$\tilde{x}_{k+1}=\left[\begin{array}[]{cc}-I_{2\times 2}&-I_{2\times
2}\end{array}\right]{\tt ReLU}\left(\left[\begin{array}[]{cc}I_{2\times
2}&\mymathbb{0}_{2\times 2}\\\ \mymathbb{0}_{2\times 2}&I_{2\times
2}\end{array}\right]\left[\begin{array}[]{c}x_{k}\\\
u_{k}\end{array}\right]+\left[\begin{array}[]{c}-50\cdot\mymathbb{1}_{2\times
1}\\\ -50\cdot\mymathbb{1}_{2\times
1}\end{array}\right]\right)+100\cdot\mymathbb{1}_{2\times 1},\vspace{-0.5em}$

which is equivalent to $\tilde{x}_{k+1}=\tilde{f}(x_{k},u_{k})=x_{k}+u_{k}$
given $x_{k}\in\mathcal{X}$ and $u_{k}\in\mathcal{U}$. Finally, the
uncertainty bounds are set to
$\epsilon^{x}=\epsilon^{u}=\epsilon^{y}=[0.05,\;0.05]^{T}$. We develop a
reachability-guided RRT planner similar to [Shkolnik et al.(2009)Shkolnik,
Walter, and Tedrake], combined with the Dijkstra algorithm, to generate a path
of reference states $x^{r}$. Different from the classic RRT, the reachability-
guided RRT incorporates dynamics as a constraint to extend the edges of the
tree. This tree $\mathcal{T}(\\{x^{(i)}\\},\\{E^{(ij)}\\})$ comprises nodes
$x^{(i)}\in\mathcal{X}_{s}$ and directed edges $E^{(ij)}$. The edge $E^{(ij)}$
connects node $x^{(i)}$ to node $x^{(j)}$ and implies that there exists a
control $u_{k}\in\mathcal{U}$ such that
$x^{(j)}=\tilde{f}(x^{(i)},u_{k})\in\mathcal{X}_{s}$, i.e.,
$x^{(j)}\in\tilde{\mathcal{F}}(\\{x^{(i)}\\},\mathcal{U})\cap\mathcal{X}_{s}$.
As shown in Figure 2, the algorithm is initialized with an initial node
$x_{0}$ and is terminated if there exists a node $x^{(i)}\in\mathcal{T}$ such
that $x_{g}\in\tilde{\mathcal{F}}(\\{x^{(i)}\\},\mathcal{U})$. At each time
step, we take the first node $x^{(i)}$ in the shortest path as the reference
state $x^{r}$ in (17). We solve the optimization (17) and apply the resulting
control $\tilde{u}_{k}$ to the actual dynamic model $f$. Then, we remove
$x^{(i)}$ from the path if
$x^{(i)}\in[\underline{x}_{k+1},\overline{x}_{k+1}]$ and the navigation
terminates when $x_{g}\in[\underline{x}_{k+1},\overline{x}_{k+1}]$. As shown
in Figure 2, after taking control $\tilde{u}_{k}$, the resulting next states
$x_{k+1}$ in blue asterisks are within the reachable set
$[\underline{x}_{k+1},\overline{x}_{k+1}]$ in green boxes and the boxes
$[\underline{x}_{k+1},\overline{x}_{k+1}]$ are collision-free, which
empirically validates the Proposition 3.5.

##  5 Vehicle Navigation and Set-theoretical Localization

Figure 3: Schematic of navigating a vehicle through a maze with tracking
controller (17) and with a set-theoretic localization algorithm. The zoom-in
view, at the top-left, demonstrates that the set-theoretic localization
algorithm provides estimates of the vehicle body and orientation that are
guaranteed to contain the actual ones. The tracking controller (17) leverages
this information, together with an NN-learned vehicle dynamics model, to avoid
obstacles.

As shown in Figure 3, we consider a front-wheel drive vehicle of width $2\;\rm
m$. The length of the vehicle wheelbase is $l=5\;\rm m$. We adopt the vehicle
kinematics model from [Li et al.(2023)Li, Li, Li, Girard, and Kolmanovsky],
which admits the form

$x_{k+1}=f(x_{k},u_{k})=\left[\begin{array}[]{c}p_{x,k}+v_{k}dt\cos{\theta_{k}}\cos{\delta_{k}}\\\
p_{x,k}+v_{k}dt\sin{\theta_{k}}\cos{\delta_{k}}\\\
\theta_{k}+v_{k}dt/l\sin{\delta_{k}}\end{array}\right],\vspace{-0.5em}$

where $x_{k}=[p_{x,k}\;p_{y,k}\;\theta_{k}]^{T}$ is the state vector,
$(p_{x,k},\;p_{y,k})$ in meters are the coordinates of the center of the
vehicle rear wheel axis, and $\theta_{k}\in[-\pi,\pi]$ is the vehicle
orientation; $u_{k}=[v_{k}\;\delta_{k}]^{T}$ is the control vector, $v_{k}$ in
$\rm m/s$ is the vehicle longitudinal speed, and $\delta_{k}$ in $\rm rad$ is
the vehicle steering angle; $dt=0.1\;\rm sec$ is the sampling period. We use
the sets
$\mathcal{U}=\left\\{u_{k}:v_{k}\in[2,5],\;\delta_{k}\in[-0.6,0.6]\right\\}$,
$\mathcal{X}=\left\\{x_{k}:p_{x,k},p_{y,k}\in[-50,50],\;\theta_{k}\in[-\pi,\pi]\right\\}$
together with obstacles $\mathcal{X}_{u}^{(i)}$ visualized in grey boxes in
Figure 3. For the uncertainties, we assume the actuator disturbance
$\epsilon^{u}=[0.01\;\rm m/s,\;0.5\;\rm deg]^{T}$. We densely sample a dataset
$\mathcal{D}=\\{(x^{(i)}_{k},u^{(i)}_{k},x^{(i)}_{k+1})\\}_{i}$ from
$\mathcal{X}\times\mathcal{U}$ for NN training and quantifying the NN
prediction error. Based on Pytorch library([Paszke et al.(2019)Paszke, Gross,
Massa, Lerer, Bradbury, Chanan, Killeen, Lin, Gimelshein, Antiga, et al.]), we
train an NN $\tilde{f}$ using the Stochastic Gradient Descent algorithm and
dataset $\mathcal{D}$ to minimize the mean-squared error $\left\lVert
x^{(i)}_{k+1}-\tilde{f}(x^{(i)}_{k},u^{(i)}_{k})\right\rVert^{2}_{2}$, and the
prediction error is equal to $\epsilon^{x}=[0.02\;\rm m,\;0.02\;\rm
m,\;1.5\;\rm deg]^{T}$. We quantify the prediction error as the maximum value
of the empirical absolute error, i.e.,

$\epsilon^{x}=\max\left\\{\epsilon\in\mathbb{R}^{3}:\epsilon=\left|x^{(i)}_{k+1}-\tilde{f}(x^{(i)}_{k},u^{(i)}_{k})\right|,\;(x^{(i)}_{k},u^{(i)}_{k},x^{(i)}_{k+1})\in\mathcal{D}\right\\},\vspace{-0.5em}$

and the function $\max$ is applied element-wise. The theoretical properties of
uncertainty-bound quantification from samples are discussed in [Dean et
al.(2020)Dean, Matni, Recht, and Ye] and are beyond the scope of this work. At
each time step $k$, we apply the set-theoretic localization algorithm
presented in [Li et al.(2023)Li, Li, Li, Girard, and Kolmanovsky] to generate
an uncertainty polygon $P_{xy}$ that contains the actual vehicle position
$[p_{x,k},\;p_{y,k}]^{T}$ and an uncertainty interval $P_{\theta}$ that
contains the actual vehicle orientation $\theta_{k}$, i.e.,
$[p_{x,k},\;p_{y,k}]^{T}\in P_{xy}$ and $\theta_{k}\in P_{\theta}$.
Subsequently, we can derive the smallest hypercube $P$ that over-bounds
$P_{xy}\times P_{\theta}$, and the measurement error $\epsilon^{y}$ is
quantified as half of the sizes of $P$. The algorithm can also produce a
polytope estimation of the vehicle body as demonstrated in Figure 3. Then,
akin to the process in Section 4, we solve the optimization problem (17) and
use the resulting control to navigate the vehicle through the maze. As shown
in Figure 3, we can also observe that the vehicle is collision-free, which is
guaranteed formally by Proposition 3.5.

##  6 Conclusion and Future Work

In this paper, we developed an approach for robust reference tracking that
leveraged a learned NN model to control the actual dynamics. We considered
both bounded intrinsic and extrinsic uncertainties from the controller and
other system modules, respectively. We transcripted the resulting variable
uncertainty set propagation through NN using a MILP. We provided formal proof
that the proposed MILP can render the overall system safe considering all
possible actuator disturbances, measurement noise, and prediction errors
within their corresponding bounded sets. We test the proposed method in
navigation and obstacle avoidance of omnidirectional robots and vehicles in
simulations. Future work will include the investigation of tightening the
over-estimation by the hypercubes, and an extension to ensure the recursive
feasibility and stability of the proposed MILP.

## References

  * [Aswani et al.(2013)Aswani, Gonzalez, Sastry, and Tomlin] Anil Aswani, Humberto Gonzalez, S Shankar Sastry, and Claire Tomlin.  Provably safe and robust learning-based model predictive control.  _Automatica_ , 49(5):1216–1226, 2013.
  * [Barnes et al.(2017)Barnes, Maddern, and Posner] Dan Barnes, Will Maddern, and Ingmar Posner.  Find your own way: Weakly-supervised segmentation of path proposals for urban autonomy.  In _2017 IEEE International Conference on Robotics and Automation (ICRA)_ , pages 203–210. IEEE, 2017.
  * [Berkenkamp et al.(2016)Berkenkamp, Moriconi, Schoellig, and Krause] Felix Berkenkamp, Riccardo Moriconi, Angela P Schoellig, and Andreas Krause.  Safe learning of regions of attraction for uncertain, nonlinear systems with gaussian processes.  In _Conference on Decision and Control_ , pages 4661–4666. IEEE, 2016.
  * [Dai et al.(2017)Dai, Gan, Xia, and Zhan] Liyun Dai, Ting Gan, Bican Xia, and Naijun Zhan.  Barrier certificates revisited.  _Journal of Symbolic Computation_ , 80:62–86, 2017.
  * [Dean and Recht(2021)] Sarah Dean and Benjamin Recht.  Certainty equivalent perception-based control.  In _Learning for Dynamics and Control_ , pages 399–411. PMLR, 2021.
  * [Dean et al.(2020)Dean, Matni, Recht, and Ye] Sarah Dean, Nikolai Matni, Benjamin Recht, and Vickie Ye.  Robust guarantees for perception-based control.  In _Learning for Dynamics and Control_ , pages 350–360. PMLR, 2020.
  * [Dean et al.(2021)Dean, Taylor, Cosner, Recht, and Ames] Sarah Dean, Andrew Taylor, Ryan Cosner, Benjamin Recht, and Aaron Ames.  Guaranteeing safety of learned perception modules via measurement-robust control barrier functions.  In _Conference on Robot Learning_ , pages 654–670. PMLR, 2021.
  * [Gal and Ghahramani(2016)] Yarin Gal and Zoubin Ghahramani.  Dropout as a bayesian approximation: Representing model uncertainty in deep learning.  In _international conference on machine learning_ , pages 1050–1059. PMLR, 2016.
  * [Hewing et al.(2020)Hewing, Wabersich, Menner, and Zeilinger] Lukas Hewing, Kim P Wabersich, Marcel Menner, and Melanie N Zeilinger.  Learning-based model predictive control: Toward safe learning in control.  _Annual Review of Control, Robotics, and Autonomous Systems_ , 3:269–296, 2020.
  * [Hu et al.(2020)Hu, Fazlyab, Morari, and Pappas] Haimin Hu, Mahyar Fazlyab, Manfred Morari, and George J Pappas.  Reach-sdp: Reachability analysis of closed-loop systems with neural network controllers via semidefinite programming.  In _2020 59th IEEE Conference on Decision and Control (CDC)_ , pages 5929–5934. IEEE, 2020.
  * [Hu et al.(2023)Hu, Yang, Chen, Li, Sima, Zhu, Chai, Du, Lin, Wang, et al.] Yihan Hu, Jiazhi Yang, Li Chen, Keyu Li, Chonghao Sima, Xizhou Zhu, Siqi Chai, Senyao Du, Tianwei Lin, Wenhai Wang, et al.  Planning-oriented autonomous driving.  In _Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition_ , pages 17853–17862, 2023.
  * [LeCun et al.(2015)LeCun, Bengio, and Hinton] Yann LeCun, Yoshua Bengio, and Geoffrey Hinton.  Deep learning.  _nature_ , 521(7553):436–444, 2015.
  * [Li et al.(2021)Li, Du, Zeng, and Jenkins] Xiao Li, Yidong Du, Zhen Zeng, and Odest Chadwicke Jenkins.  Seannet: Semantic understanding network for localization under object dynamics.  _arXiv preprint arXiv:2110.02276_ , 2021.
  * [Li et al.(2023)Li, Li, Li, Girard, and Kolmanovsky] Xiao Li, Yutong Li, Nan Li, Anouck Girard, and Ilya Kolmanovsky.  Set-theoretic localization for mobile robots with infrastructure-based sensing.  _Advanced Control for Applications: Engineering and Industrial Systems_ , 5(1):e117, 2023.
  * [Lofberg(2004)] Johan Lofberg.  Yalmip: A toolbox for modeling and optimization in matlab.  In _2004 IEEE international conference on robotics and automation (IEEE Cat. No. 04CH37508)_ , pages 284–289. IEEE, 2004.
  * [Paszke et al.(2019)Paszke, Gross, Massa, Lerer, Bradbury, Chanan, Killeen, Lin, Gimelshein, Antiga, et al.] Adam Paszke, Sam Gross, Francisco Massa, Adam Lerer, James Bradbury, Gregory Chanan, Trevor Killeen, Zeming Lin, Natalia Gimelshein, Luca Antiga, et al.  Pytorch: An imperative style, high-performance deep learning library.  _Advances in Neural Information Processing Systems_ , 32, 2019.
  * [Roddick and Cipolla(2020)] Thomas Roddick and Roberto Cipolla.  Predicting semantic map representations from images using pyramid occupancy networks.  In _Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition_ , pages 11138–11147, 2020.
  * [Shkolnik et al.(2009)Shkolnik, Walter, and Tedrake] Alexander Shkolnik, Matthew Walter, and Russ Tedrake.  Reachability-guided sampling for planning under differential constraints.  In _2009 IEEE International Conference on Robotics and Automation_ , pages 2859–2865. IEEE, 2009.
  * [Tjeng et al.(2017)Tjeng, Xiao, and Tedrake] Vincent Tjeng, Kai Xiao, and Russ Tedrake.  Evaluating robustness of neural networks with mixed integer programming.  _arXiv preprint arXiv:1711.07356_ , 2017.
  * [Tong et al.(2023)Tong, Dawson, and Fan] Mukun Tong, Charles Dawson, and Chuchu Fan.  Enforcing safety for vision-based controllers via control barrier functions and neural radiance fields.  In _2023 IEEE International Conference on Robotics and Automation (ICRA)_ , pages 10511–10517. IEEE, 2023.
  * [Wang et al.(2018)Wang, Han, and Egerstedt] Li Wang, Dongkun Han, and Magnus Egerstedt.  Permissive barrier certificates for safe stabilization using sum-of-squares.  In _2018 Annual American Control Conference (ACC)_ , pages 585–590. IEEE, 2018.
  * [Wei and Liu(2021)] Tianhao Wei and Changliu Liu.  Safe control with neural network dynamic models.  _arXiv preprint arXiv:2110.01110_ , 2021.

##  Appendix A Formal Safety Guarantee

The proofs of Proposition 3.1, 3.2, 3.3 are presented in A.1, A.2, A.3,
respectively. Subsequently, the hypercube defined by the lower and upper
bounds $\underline{x}_{k+1},\overline{x}_{k+1}$, as a result of the MILP
problem (17), contains the actual state $x_{k+1}$, and the hypercube is within
the safe set $\mathcal{X}_{s}$, i.e.,
$x_{k+1}\in[\underline{x}_{k+1},\overline{x}_{k+1}]\subset\mathcal{X}_{s}$.
This proves Proposition 3.5.

###  A.1 Input Constraints

######  Proof A.1.

Based on assumption in (3) and $x_{k}\in\mathcal{X}_{s}\subset\mathcal{X}$,
it’s obvious that
$x_{k}\in\mathcal{X}_{k}=(y_{k}\ominus\mathcal{W}_{y})\cap\mathcal{X}_{s}$ by
the definition of Pontryagin difference, thereby we have
$x_{k}\in\mathcal{X}_{k}\subset(y_{k}\ominus\mathcal{W}_{y})\cap\mathcal{X}$.
As defined in (9) together with the hypercube assumptions in (6b), (8a), the
following statement can be proved

$\forall
x\in(y_{k}\ominus\mathcal{W}_{y})\cap\mathcal{X},\;a_{0,1:n_{x}}=\max\\{\underline{x},y_{k}-\epsilon^{y}\\}\leq
x\leq\min\\{\overline{x},y_{k}+\epsilon^{y}\\}=b_{0,1:n_{x}}.$

Therefore, we can conclude that

$x_{k}\in\left[a_{0,1:n_{x}},\;b_{0,1:n_{x}}\right].$ (18)

Similarly, based on (4) and $u_{k}\in\mathcal{U}$, we can show that
$u_{k}\in\mathcal{U}_{k}$ and
$\max\\{\underline{u},\tilde{u}_{k}-\epsilon^{u}\\}\leq
u\leq\min\\{\overline{u},\tilde{u}_{k}+\epsilon^{u}\\}$ for all
$u\in\mathcal{U}_{k}$. Next, we need to show that the constraints (10) imply
$a_{0,(n_{x}+1):(n_{x}+n_{u})}=\max\\{\underline{u},\tilde{u}_{k}-\epsilon^{u}\\}$
and
$b_{0,(n_{x}+1):(n_{x}+n_{u})}=\min\\{\overline{u},\tilde{u}_{k}+\epsilon^{u}\\}$
from which we can conclude

$u_{k}\in\left[a_{0,(n_{x}+1):(n_{x}+n_{u})},\;b_{0,(n_{x}+1):(n_{x}+n_{u})}\right].$
(19)

If $\delta^{a}_{j}=1$, the constraints (10) yield the following inequalities

$\left\\{\begin{array}[]{l}a_{0,j}\geq\underline{u}_{j}\\\
a_{0,j}\geq\tilde{u}_{k,j}-\epsilon^{u}_{j}\\\ a_{0,j}\leq\underline{u}_{j}\\\
a_{0,j}\leq\tilde{u}_{k,j}-\epsilon^{u}_{j}+\max\\{\epsilon^{u}_{j},\overline{u}_{j}-\underline{u}_{j}-\epsilon^{u}_{j}\\}\end{array}\right.$

where the first three inequalities are equivalent to
$a_{0,j}=\underline{u}_{j},\;\underline{u}_{j}\geq\tilde{u}_{k,j}-\epsilon^{u}_{j}$.
In this case, we have
$a_{0,j}=\max\\{\underline{u}_{j},\tilde{u}_{k,j}-\epsilon^{u}_{j}\\}$ and the
fourth inequality is valid since
$\underline{u}_{j}\leq\tilde{u}_{k,j}\leq\overline{u}_{j}$ and
$a_{0,j}-\tilde{u}_{k,j}+\epsilon^{u}_{j}\leq\underline{u}_{j}-\underline{u}_{j}+\epsilon^{u}_{j}\leq\epsilon^{u}_{j}\leq\max\\{\epsilon^{u}_{j},\overline{u}_{j}-\underline{u}_{j}-\epsilon^{u}_{j}\\}$.

If $\delta^{a}_{j}=0$, the constraints (10) produce the following inequalities

$\left\\{\begin{array}[]{l}a_{0,j}\geq\underline{u}_{j}\\\
a_{0,j}\geq\tilde{u}_{k,j}-\epsilon^{u}_{j}\\\
a_{0,j}\leq\tilde{u}_{k,j}-\epsilon^{u}_{j}\\\
a_{0,j}\leq\underline{u}_{j}+\max\\{\epsilon^{u}_{j},\overline{u}_{j}-\underline{u}_{j}-\epsilon^{u}_{j}\\}\end{array}\right.$

where the first three inequalities are equivalent to
$a_{0,j}=\tilde{u}_{k,j}-\epsilon^{u}_{j},\;\tilde{u}_{k,j}-\epsilon^{u}_{j}\geq\underline{u}_{j}$.
Again, we can show
$a_{0,j}=\max\\{\underline{u}_{j},\tilde{u}_{k,j}-\epsilon^{u}_{j}\\}$ and the
fourth inequality is also feasible since
$a_{0,j}-\underline{u}_{j}=\tilde{u}_{k,j}-\epsilon^{u}_{j}-\underline{u}_{j}\leq\overline{u}_{j}-\epsilon^{u}_{j}-\underline{u}_{j}\leq\max\\{\epsilon^{u}_{j},\overline{u}_{j}-\underline{u}_{j}-\epsilon^{u}_{j}\\}$.

Hence, we can conclude that
$a_{0,j}=\max\\{\underline{u}_{j},\tilde{u}_{k,j}-\epsilon^{u}_{j}\\}$ for
$j=(n_{x}+1),\dots,(n_{x}+n_{u})$ that is equivalent to
$a_{0,(n_{x}+1):(n_{x}+n_{u})}=\max\\{\underline{u},\tilde{u}_{k}-\epsilon^{u}\\}$
and the proof of
$b_{0,(n_{x}+1):(n_{x}+n_{u})}=\min\\{\overline{u},\tilde{u}_{k}+\epsilon^{u}\\}$
resembles the discussion above. Eventually, based on (18) and (19), the input
$z_{0}$ to the NN satisfies

$z_{0}=[x_{k}^{T}\;u_{k}^{T}]^{T}\in\left[a_{0},\;b_{0}\right].$ (20)

###  A.2 NN Structural Constraints

######  Proof A.2.

We first provide proof of the following statement

$\hat{z}_{i}=W^{(i)}z_{i-1}+b^{(i)}\in\left[\hat{a}_{i},\;\hat{b}_{i}\right],\;z_{i}=\max\\{0,\;\hat{z}_{i}\\}\in\left[a_{i},\;b_{i}\right],\;\text{given
}z_{i-1}\in\left[a_{i-1},\;b_{i-1}\right].$

Then, from the results $z_{0}\in\left[a_{0},\;b_{0}\right]$ in (20), we can
inductively prove the same argument for $i=1,\dots,\ell-1$.

Given $z_{i-1}\in\left[a_{i-1},\;b_{i-1}\right]$, the $j^{th}$ element in
$\hat{z}_{i}$ is $\hat{z}_{i,j}=w^{(i)}_{j}z_{i-1}+b^{(i)}_{j}$ where
$w^{(i)}_{j}=[w_{1}\dots w_{q}\dots w_{n_{i-1}}]$ and satisfies the following
inequalities

$\begin{array}[]{c}\sum\limits_{q=1}^{n_{i-1}}\left(\llbracket w_{q}\geq
0\rrbracket\cdot\left(w_{q}a_{i-1,q}\right)+\llbracket
w_{q}<0\rrbracket\cdot\left(w_{q}b_{i-1,q}\right)\right)\\\
\leq\hat{z}_{i,j}-b^{(i)}_{j}\leq\\\
\sum\limits_{q=1}^{n_{i-1}}\left(\llbracket w_{q}\geq
0\rrbracket\cdot\left(w_{q}b_{i-1,q}\right)+\llbracket
w_{q}<0\rrbracket\cdot\left(w_{q}a_{i-1,q}\right)\right),\end{array}$

where the Iverson bracket $\llbracket\cdot\rrbracket$ takes the value of 1 if
the statement inside is true and 0 otherwise. We can derive the following
equalities from the constraints in (11)

$\displaystyle\hat{a}_{i,j}-b^{(i)}_{j}$
$\displaystyle=\sum\limits_{q=1}^{n_{i-1}}\left(\llbracket w_{q}\geq
0\rrbracket\cdot\left(w_{q}a_{i-1,q}\right)+\llbracket
w_{q}<0\rrbracket\cdot\left(w_{q}b_{i-1,q}\right)\right),$
$\displaystyle\hat{b}_{i,j}-b^{(i)}_{j}$
$\displaystyle=\sum\limits_{q=1}^{n_{i-1}}\left(\llbracket w_{q}\geq
0\rrbracket\cdot\left(w_{q}b_{i-1,q}\right)+\llbracket
w_{q}<0\rrbracket\cdot\left(w_{q}a_{i-1,q}\right)\right),$

which implies

$\hat{z}_{i,j}\in\left[\hat{a}_{i,j},\;\hat{b}_{i,j}\right],\;j=1,\dots,n_{i}.$
(21)

Afterward, given
$\underline{\hat{z}}_{i,j}\leq\hat{a}_{i,j}\leq\hat{z}_{i,j}\leq\hat{b}_{i,j}\leq\overline{\hat{z}}_{i,j}$
in (11) and constraints in (12), we need to prove $a_{i,j}\leq z_{i,j}\leq
b_{i,j}$. If $\delta^{--}_{i,j}=1$, the constraints in (12) are reduced to

$\left\\{\begin{array}[]{l}b_{i,j}\geq a_{i,j}\geq 0\\\ \hat{a}_{i,j}\leq
a_{i,j}\leq 0\\\ \hat{b}_{i,j}\leq b_{i,j}\leq 0\\\
a_{i,j}\leq\hat{a}_{i,j}-\underline{\hat{z}}_{i,j}\\\
b_{i,j}\leq\hat{b}_{i,j}-\underline{\hat{z}}_{i,j}\end{array}\right..$

The first three inequalities are equivalent to $a_{i,j}=b_{i,j}=0$,
$\hat{a}_{i,j}\leq 0$, $\hat{b}_{i,j}\leq 0$ which induce
$\hat{z}_{i,j}\leq\hat{b}_{i,j}\leq 0$ and
$z_{i,j}=\max\\{0,\;\hat{z}_{i,j}\\}=0$, thereby we have $0=a_{i,j}\leq
z_{i,j}\leq b_{i,j}=0$. Meanwhile, the last two inequalities hold since
$\underline{\hat{z}}_{i,j}\leq\hat{a}_{i,j}\leq\hat{z}_{i,j}\leq\hat{b}_{i,j}\leq\overline{\hat{z}}_{i,j}$.

If $\delta^{-+}_{i,j}=1$, the constraints in (12) are reduced to

$\left\\{\begin{array}[]{l}b_{i,j}\geq a_{i,j}\geq 0\\\ \hat{a}_{i,j}\leq
a_{i,j}\leq 0\\\ \hat{b}_{i,j}\leq b_{i,j}\leq\hat{b}_{i,j}\\\
a_{i,j}\leq\hat{a}_{i,j}-\underline{\hat{z}}_{i,j}\\\
b_{i,j}\leq\overline{\hat{z}}_{i,j}\end{array}\right..$

The first three inequalities imply that $a_{i,j}=0$, $\hat{a}_{i,j}\leq 0$ and
$b_{i,j}=\hat{b}_{i,j}$, $\hat{b}_{i,j}\geq 0$ from which we can show that
$0\leq z_{i,j}\leq\hat{b}_{i,j}$. Thus, we have $a_{i,j}=0\leq
z_{i,j}\leq\hat{b}_{i,j}=b_{i,j}$ and the last two inequalities are also
feasible.

If $\delta^{++}_{i,j}=1$, the constraints in (12) are reduced to

$\left\\{\begin{array}[]{l}b_{i,j}\geq a_{i,j}\geq 0\\\ \hat{a}_{i,j}\leq
a_{i,j}\leq\hat{a}_{i,j}\\\ \hat{b}_{i,j}\leq b_{i,j}\leq\hat{b}_{i,j}\\\
a_{i,j}\leq\overline{\hat{z}}_{i,j}\\\
b_{i,j}\leq\overline{\hat{z}}_{i,j}\end{array}\right..$

The first three inequalities imply that $a_{i,j}=\hat{a}_{i,j}$,
$\hat{a}_{i,j}\geq 0$ and $b_{i,j}=\hat{b}_{i,j}$, $\hat{b}_{i,j}\geq 0$ from
which we can show $\hat{a}_{i,j}\leq z_{i,j}\leq\hat{b}_{i,j}$. Thus, we have
$a_{i,j}=\hat{a}_{i,j}\leq z_{i,j}\leq\hat{b}_{i,j}=b_{i,j}$ and the last two
inequalities are valid. To this end, we have proven that

$z_{i,j}\in\left[a_{i,j},\;b_{i,j}\right],\;j=1,\dots,n_{i}.$ (22)

Indeed, as discussed at the beginning of this section, we can show that

$\hat{z}_{i}\in\left[\hat{a}_{i},\;\hat{b}_{i}\right],\;z_{i}\in\left[a_{i},\;b_{i}\right],\;i=1,\dots,\ell-1$
(23)

with (20). Identical to the proof of (21), given
$z_{\ell-1}\in\left[a_{\ell-1},\;b_{\ell-1}\right]$ in (23), constraints in
(13) and $\tilde{x}_{k+1}=W^{(\ell)}z_{\ell-1}+b^{(\ell)}$ in (1), we can
demonstrate that

$\tilde{x}_{k+1}\in\tilde{\mathcal{F}}(\mathcal{X}_{k},\mathcal{U}_{k})\subseteq\left[a_{k+1},\;b_{k+1}\right].$
(24)

###  A.3 Output Constraints

######  Proof A.3.

With the assumptions in (2) and (6a), we can show
$x_{k+1}\in\left[a_{k+1},\;b_{k+1}\right]\oplus\mathcal{W}_{x}$ which can be
further derived to

$x_{k+1}\in\left[\underline{x}_{k+1},\;\overline{x}_{k+1}\right]\subset\mathcal{X},$
(25)

based on the constraints in (14). Meanwhile, the constraints in (15) are
equivalent to the following statement

$\exists j\in\\{1,\dots,n_{x}\\},\;s.t.\;\delta^{u}_{1,j}=1\text{ or
}\delta^{u}_{2,j}=1.$

In the sequel, we show that
$\left[\underline{x}_{k+1},\;\overline{x}_{k+1}\right]\cap\mathcal{X}_{u}=\varnothing$
with constraints in (15) if $\delta^{u}_{1,j}=1$ which the proof resembles in
case of $\delta^{u}_{2,j}=1$. If $\delta^{u}_{1,j}=1$, the constraints in (15)
can be reformulated into

$\underline{x}_{j}\leq\underline{x}_{k+1,j}\leq\overline{x}_{k+1,j}\leq\underline{x}_{u,j}.$

Suppose
$\left[\underline{x}_{k+1},\;\overline{x}_{k+1}\right]\cap\mathcal{X}_{u}\neq\varnothing$,
therefore there exists $x^{*}=[x^{*}_{1}\cdots x^{*}_{j}\cdots
x^{*}_{n_{x}}]^{T}\in\mathbb{R}^{n_{x}}$ such that
$x^{*}\in\left[\underline{x}_{k+1},\;\overline{x}_{k+1}\right]\cap\mathcal{X}_{u}$.
Then, for some index $j$ and $1\leq j\leq n_{x}$, we have
$\underline{x}_{u,j}\leq x^{*}_{j}\leq\overline{x}_{u,j}$ and
$\underline{x}_{k+1,j}\leq x^{*}_{j}\leq\overline{x}_{k+1,j}$ which violates
$\overline{x}_{k+1,j}\leq\underline{x}_{u,j}$. By contradiction, the following
is true

$\left[\underline{x}_{k+1},\;\overline{x}_{k+1}\right]\cap\mathcal{X}_{u}=\varnothing.$
(26)

With equations (25) and (26), we demonstrate that
$x_{k+1}\in\mathcal{F}(\mathcal{X}_{k},\mathcal{U}_{k})\subseteq[\underline{x}_{k+1},\overline{x}_{k+1}]\subset\mathcal{X}_{s}$
which complete the proof for Proposition 3.3.