#  Boundary and Taxonomy of Integrator Reach Sets

Shadi Haddad, Abhishek Halder Shadi Haddad and Abhishek Halder are with the
Department of Applied Mathematics, University of California, Santa Cruz, CA
95064, USA<EMAIL_ADDRESS>
###### Abstract

Over-approximating the forward reach sets of controlled dynamical systems
subject to set-valued uncertainties is a common practice in systems-control
engineering for the purpose of performance verification. However, specific
algebraic and topological results for the geometry of such sets are rather
uncommon even for simple linear systems such as the integrators. This work
explores the geometry of the forward reach set of the integrator dynamics
subject to box-valued uncertainties in its control inputs. Our contribution
includes derivation of a closed-form formula for the support functions of
these sets. This result, then enables us to deduce the parametric as well as
the implicit equations describing the exact boundaries of these reach sets.
Specifically, the implicit equations for the bounding hypersurfaces are shown
to be given by vanishing of certain Hankel determinants. Finally, it is
established that these sets are semialgebraic as well as translated zonoids.
Such results should be useful to benchmark existing reach set over-
approximation algorithms, and to help design new algorithms for the same.

##  I Introduction

Consider the forward reach set $\mathcal{R}$ for the integrator dynamics in
Brunovsky normal form with $d$ states and $m$ inputs given by

$\displaystyle\dot{\bm{x}}=\underbrace{\mathrm{blkdiag}\left(\bm{A}_{1},\ldots,\bm{A}_{m}\right)}_{=:\bm{A}}\bm{x}+\underbrace{\mathrm{blkdiag}\left(\bm{b}_{1},\ldots,\bm{b}_{m}\right)}_{=:\bm{B}}\bm{u},$
(1)

with a relative degree vector
$\bm{r}=(r_{1},\ldots,r_{m})^{\top}\in\mathbb{N}^{m}$ (i.e., an $m\times 1$
vector of natural numbers) where $r_{1}+\ldots+r_{m}=d$. In (1), the symbol
$\mathrm{blkdiag}\left(\cdot\right)$ denotes a block diagonal matrix whose
arguments constitute its diagonal blocks, and

$\displaystyle\bm{A}_{j}:=\left(\mathbf{0}_{r_{j}\times
1}\left|\bm{e}_{1}^{r_{j}}\right|\bm{e}_{2}^{r_{j}}|\ldots|\bm{e}_{r_{j}-1}^{r_{j}}\right),\quad\bm{b}_{j}:=\bm{e}_{r_{j}}^{r_{j}},$
(2)

for $j=1,\ldots,m$. In (2), the notation $\mathbf{0}_{r_{j}\times 1}$ stands
for the column vector of zeros of size $r_{j}\times 1$. For $k\leq\ell$, we
use $\bm{e}_{k}^{\ell}$ to denote the $k$th standard basis (column) vector in
$\mathbb{R}^{\ell}$.

As a concrete example, consider $d=5$, $m=2$, and $\bm{r}=(2,3)^{\top}$. Then
(1) takes the form

$\displaystyle\dot{\bm{x}}=\left(\begin{array}[]{@{}cc|ccc@{}}0&1&0&0&0\\\
0&0&0&0&0\\\ \hline\cr 0&0&0&1&0\\\ 0&0&0&0&1\\\ 0&0&0&0&0\\\
\end{array}\right)\bm{x}+\left(\begin{array}[]{@{}c|c@{}}0&0\\\ 1&0\\\
\hline\cr 0&0\\\ 0&0\\\ 0&1\end{array}\right)\bm{u}.$ (13)

We suppose that the initial condition $\bm{x}_{0}\in\mathbb{R}^{d}$ for (1) is
fixed, and the control input $\bm{u}\in\mathcal{U}\subset\mathbb{R}^{m}$ with
box-valued input set $\mathcal{U}$. The set-valued uncertainties in the input
$\bm{u}$ entails uncertainties in the state evolution for (1).

Specifically, the state vector
$\bm{x}(t)\in\mathcal{R}\left(\\{\bm{x}_{0}\\},t\right)\subset\mathbb{R}^{d}$,
where the set $\mathcal{R}\left(\\{\bm{x}_{0}\\},t\right)$ is the forward
reach set at time $t$ resulting from the singleton initial condition set
$\\{\bm{x}_{0}\\}$ subject to (1), i.e., the set of states the controlled
dynamics (1) may reach at time $t$ starting from the initial condition
$\bm{x}_{0}$.

Our motivation behind studying the geometry of integrator reach sets is
twofold.

  * (i)

Integrators are simple linear time invariant systems that appear frequently in
the literature (see e.g., [1, 2], [3, Ch. 3 and 4]) to benchmark the
performance of reach set over-approximation algorithms. In this context, the
lack of knowledge of exact geometry of these sets implies the lack of ground
truth against which the conservatism of different over-approximation
algorithms can be quantified. Consequently, one has to contend with graphical
or statistical performance results. The present work remedies this by deriving
the _exact_ geometry of integrator reach sets.

  * (ii)

Integrator dynamics in Brunovsky normal form also arise from differentially
flat nonlinear systems [4, 5, 6], which in turn appear frequently in vehicular
dynamics and control applications. Addressing safety and collision avoidance
for such systems in the presence of set-valued uncertainties, then requires
computing their (typically nonconvex) reach sets. Having an analytical handle
of the corresponding integrator reach sets in the normal coordinates, together
with the knowledge of coordinate transforms, could enable the design of novel
reach set computation algorithms for this class of nonlinear systems.

It is known [7] that the set $\mathcal{R}\left(\\{\bm{x}_{0}\\},t\right)$ is
compact and convex. Our prior work [8] investigated the convex geometry of the
integrator reach set for _single input case_ under the assumption that the set
$\mathcal{U}$ is a _symmetric interval_ of the form
$[-\mu,\mu]\subset\mathbb{R}$. The development herein generalizes that work–we
now allow multiple generic interval-valued inputs, i.e., $\mathcal{U}$ is a
box of the form

$\displaystyle\left[\alpha_{1},\beta_{1}\right]\times\left[\alpha_{2},\beta_{2}\right]\times\ldots\left[\alpha_{m},\beta_{m}\right]\subset\mathbb{R}^{m}.$
(14)

The box (14) can be seen as an over-approximation of arbitrary compact
$\mathcal{U}\subset\mathbb{R}^{m}$ with
$\alpha_{j}:=\underset{\bm{u}\in\mathcal{U}}{\min}\;u_{j},\quad\beta_{j}:=\underset{\bm{u}\in\mathcal{U}}{{\max}\;u_{j}}$.
The reach set resulting from (14) will be an over-approximation of that
resulting from such compact $\mathcal{U}$.

With this more generic set up, we obtain exact descriptions of the boundary
$\partial\mathcal{R}\left(\\{\bm{x}_{0}\\},t\right)$, which are novel results.
Furthermore, we show that the convex set $\mathcal{R}$ is semialgebraic,
translated zonoid, and not a spectrahedron.

One may describe $\mathcal{R}\left(\\{\bm{x}_{0}\\},t\right)$ via set-valued
operation as

$\displaystyle\mathcal{R}\left(\\{\bm{x}_{0}\\},t\right)=\exp\left(t\bm{A}\right)\bm{x}_{0}\dotplus\int_{0}^{t}\exp\left(s\bm{A}\right)\bm{B}\mathcal{U}\>{\rm{d}}s,$
(15)

with each diagonal block being upper triangular, given by

$\displaystyle\exp(\bm{A}_{j}(t-s)):=\begin{cases}\dfrac{(t-s)^{\ell-k}}{(\ell-k)!}&\text{for}\;k\leq\ell,\\\
0&\text{otherwise}.\end{cases}$ (16)

for $k,\ell=1,\ldots,r_{j}$, for each $j=1,\ldots,m$. The second summand in
(15) is a set-valued Aumann integral and $\dotplus$ stands for the Minkowski
sum. In other words, $\mathcal{R}\left(\\{\bm{x}_{0}\\},t\right)$ is simply a
translation of the set given by the second summand in (15). From here, it is
not clear if a more explicit description of $\mathcal{R}$ is possible or not.
Moreover, a description at the level of (15) is not very helpful for
computation. Even at the analysis level, it is not immediate exactly what kind
of convex set the second summand in (15) is. In this paper, we show that these
issues can be addressed in surprisingly explicit manner.

In Sec. II, we deduce the support function of the forward reach set of (1)
subject to box-valued uncertainties in its control input. The development in
Sec. II allows us to derive exact parametric (Sec. III-A) as well as implicit
equations (Sec. III-B) for the boundary of this forward reach set. In Sec. IV,
we use the results from Sec. III to show that this reach set is semialgebraic
but not a spectrahedron. Sec. V concludes the paper.

##  II Support Function of the Integrator Reach Set

The support function $h_{\mathcal{K}}:\mathbb{R}^{d}\mapsto\mathbb{R}$ of a
compact convex set $\mathcal{K}\subset\mathbb{R}^{d}$, is given by

$\displaystyle
h_{\mathcal{K}}(\bm{y}):=\underset{\bm{x}\in\mathcal{K}}{\sup}\>\\{\langle\bm{y},\bm{x}\rangle\mid\bm{y}\in\mathbb{R}^{d}\\},$
(17)

where $\langle\cdot,\cdot\rangle$ denotes the standard Euclidean inner
product. Geometrically, $h_{\mathcal{K}}(\bm{y})$ gives the signed distance of
the supporting hyperplane of $\mathcal{K}$ with outer normal vector $\bm{y}$,
measured from the origin. the distance is negative if an only if $\bm{y}$
points into the open halfspace containing the origin. Furthermore, the
supporting hyperplane at $\bm{x}^{\text{bdy}}\in\partial\mathcal{K}$ is
$\langle\bm{y},\bm{x}^{\text{bdy}}\rangle=h_{\mathcal{K}}(\bm{y})$, and we can
write

$\mathcal{K}=\big{\\{}\bm{x}\in\mathbb{R}^{d}\mid\langle\bm{y},\bm{x}\rangle\leq
h_{\mathcal{K}}(\bm{y})\;\text{for all}\;\bm{y}\in\mathbb{R}^{d}\big{\\}}.$

Given a function $f:\mathbb{R}^{d}\mapsto\mathbb{R}\cup\\{+\infty\\}$, its
Legendre-Fenchel conjugate is

$\displaystyle
f^{*}(\bm{y}):=\underset{\bm{x}\in\operatorname{domain}(f)}{\sup}\>\\{\langle\bm{y},\bm{x}\rangle-f(\bm{x})\mid\bm{y}\in\mathbb{R}^{d}\\}.$
(18)

From (17)-(18), it follows that $h_{\mathcal{K}}(\bm{y})$ is the Legendre-
Fenchel conjugate of the indicator function

$\mathbf{1}_{\mathcal{K}}(\bm{x}):=\begin{cases}0&\text{if}\quad\bm{x}\in\mathcal{K},\\\
+\infty&\text{otherwise.}\end{cases}$

Thus, the support function $h_{\mathcal{K}}(\bm{y})$ uniquely determines the
set $\mathcal{K}$. Since the indicator function of a convex set is a convex
function, the biconjugate

$\displaystyle\mathbf{1}_{\mathcal{K}}^{**}(\cdot)=h_{\mathcal{K}}^{*}(\cdot)=\mathbf{1}_{\mathcal{K}}(\cdot).$
(19)

In (1), since the system matrices are block diagonal, for the input set (14),
each of the $m$ single input integrator dynamics with $r_{j}$ dimensional
state subvectors for $j=1,\ldots,m$, are decoupled from each other. Hence
$\mathcal{R}\left(\\{\bm{x}_{0}\\},t\right)\subset\mathbb{R}^{d}$ is the
Cartesian product of these single integrator reach sets:
$\mathcal{R}_{j}\left(\\{\bm{x}_{0}\\},t\right)\subset\mathbb{R}^{r_{j}}$ for
$j=1,\ldots,m$, i.e.,

$\displaystyle\mathcal{R}=\mathcal{R}_{1}\times\mathcal{R}_{2}\times\ldots\times\mathcal{R}_{m}.$
(20)

Notice that we can also write (20) as a Minkowski sum
$\mathcal{R}_{1}\dotplus\ldots\dotplus\mathcal{R}_{m}$. Since the support
function of the Minkowski sum is the sum of support functions, this allows us
to write the support function of the reach set at time $t$ as

$\displaystyle
h_{\mathcal{R}\left(\\{\bm{x}_{0}\\},t\right)}(\bm{y})=\sum_{j=1}^{m}h_{\mathcal{R}_{j}\left(\\{\bm{x}_{0}\\},t\right)}(\bm{y}_{j}).$
(21)

In the following, it will be helpful to define

$\displaystyle\bm{\xi}_{j}(s):=\left(\begin{array}[]{c}s^{r_{j}-1}/\left(r_{j}-1\right)!\\\
s^{r_{j}-2}/\left(r_{j}-2\right)!\\\ \vdots\\\ s\\\ 1\end{array}\right),\quad
j=1,\ldots,m,$ (27)

and

$\bm{\zeta}_{j}(t):=\int_{0}^{t}\bm{\xi}_{j}(s){\rm{d}}s\in\mathbb{R}^{r_{j}}.$

Also, let

$\displaystyle\mu_{j}:=\frac{\beta_{j}-\alpha_{j}}{2},\quad\nu_{j}:=\frac{\beta_{j}+\alpha_{j}}{2}.$
(28)

We have the following result.

######  Theorem 1.

The support function of the forward reach set
$\mathcal{R}(\\{\bm{x}_{0}\\},t)$ for (1) at time $t$, for the input set (14),
is

$\displaystyle h_{\mathcal{R}(\\{\bm{x}_{0}\\},t)}\left(\bm{y}\right)=$
$\displaystyle\sum_{j=1}^{m}\bigg{\langle}\bm{y}_{j},\exp\left(t\bm{A}_{j}\right)\bm{x}_{j0}+\nu_{j}\bm{\zeta}_{j}(t)\bigg{\rangle}$
$\displaystyle+\mu_{j}\int_{0}^{t}|\langle\bm{y}_{j},\bm{\xi}_{j}(s)\rangle|\>{\rm{d}}s,$
(29)

where $\bm{x}_{j0}$ is the $r_{j}$ dimensional subvector of the initial
condition $\bm{x}_{0}\in\mathbb{R}^{d}$ corresponding to the relative degree
component $r_{j}$, where $j=1,\ldots,m$.

###### Proof.

Since support function is distributive over sum, we have

$\displaystyle
h_{\mathcal{R}_{j}(\\{\bm{x}_{0}\\},t)}\left(\bm{y}_{j}\right)=$
$\displaystyle\langle\bm{y}_{j},\exp\left(t\bm{A}_{j}\right)\bm{x}_{j0}\rangle$
$\displaystyle+h_{\int_{0}^{t}\bm{\xi}_{j}(s)[\alpha_{j},\beta_{j}]{\rm{d}}s}(\bm{y}_{j}).$
(30)

The second summand in the right hand side above, is the support function of
the set $\int_{0}^{t}\bm{\xi}_{j}(s)[\alpha_{j},\beta_{j}]{\rm{d}}s$, which we
simplify next.

From the definition of support function, we have

$\displaystyle
h_{\bm{\xi}_{j}(s)[\alpha_{j},\beta_{j}]}(\bm{y}_{j})=\underset{u_{j}\in[\alpha_{j},\beta_{j}]}{\sup}\>\langle\bm{y}_{j},\bm{\xi}_{j}(s)u_{j}\rangle.$
(31)

Since $\mathcal{U}$ is compact, the optimizer $u_{j}^{\text{opt}}$ in (31) is
a maximizer, and is given by

$\displaystyle
u_{j}^{\text{opt}}=\begin{cases}\beta_{j}&\text{for}\quad\langle\bm{y}_{j},\bm{\xi}_{j}(s)\rangle\geq
0,\\\
\alpha_{j}&\text{for}\quad\langle\bm{y}_{j},\bm{\xi}_{j}(s)\rangle<0.\end{cases}$
(32)

In terms of the Heaviside step function $H(\cdot)$, we can write

$\displaystyle u_{j}^{\text{opt}}$
$\displaystyle=\alpha_{j}+(\beta_{j}-\alpha_{j})H(\langle\bm{y}_{j},\bm{\xi}_{j}(s)\rangle)$
$\displaystyle=\alpha_{j}+(\beta_{j}-\alpha_{j})\times\frac{1}{2}\left(1+{\rm{sgn}}\left(\langle\bm{y}_{j},\bm{\xi}_{j}(s)\rangle\right)\right),$

where ${\rm{sgn}}(\cdot)$ denotes the signum function. Therefore,

$\displaystyle\\!\\!\\!h_{\bm{\xi}_{j}(s)[\alpha_{j},\beta_{j}]}(\bm{y}_{j})$
$\displaystyle=u_{j}^{\text{opt}}\langle\bm{y}_{j},\bm{\xi}_{j}(s)\rangle$
$\displaystyle=\nu_{j}\langle\bm{y}_{j},\bm{\xi}_{j}(s)\rangle+\mu_{j}|\langle\bm{y}_{j},\bm{\xi}_{j}(s)\rangle|,$
(33)

since $|x|=x\>{\rm{sgn}}(x)$ for any real $x$.

Applying [8, Proposition 1] to (33), we obtain

$\displaystyle
h_{\int_{0}^{t}\bm{\xi}_{j}(s)[\alpha_{j},\beta_{j}]{\rm{d}}s}(\bm{y}_{j})=\int_{0}^{t}\\!\\!h_{\bm{\xi}_{j}(s)[\alpha_{j},\beta_{j}]}(\bm{y}_{j}){\rm{d}}s$
$\displaystyle=\bigg{\langle}\bm{y}_{j},\nu_{j}\bm{\zeta}_{j}(t)\bigg{\rangle}+\mu_{j}\int_{0}^{t}|\langle\bm{y}_{j},\bm{\xi}_{j}(s)\rangle|\>{\rm{d}}s,$

which upon substituting back in (30), yields

$\displaystyle
h_{\mathcal{R}_{j}(\\{\bm{x}_{0}\\},t)}\left(\bm{y}_{j}\right)=$
$\displaystyle\bigg{\langle}\bm{y}_{j},\exp\left(t\bm{A}_{j}\right)\bm{x}_{j0}+\nu_{j}\bm{\zeta}_{j}(t)\bigg{\rangle}$
$\displaystyle+\mu_{j}\int_{0}^{t}|\langle\bm{y}_{j},\bm{\xi}_{j}(s)\rangle|\>{\rm{d}}s.$
(34)

Combining (34) with (21), we arrive at (29). ∎

######  Remark 1.

The formula (34) has a clear geometric meaning. The integral term is the
support function of a scaled zonoid, see [8, Sec. III.A.1]. The set
$\mathcal{R}_{j}(\\{\bm{x}_{0}\\},t)$ is then a translated zonoid, with a
time-varying translation of amount
$\exp\left(t\bm{A}_{j}\right)\bm{x}_{j0}+\nu_{j}\bm{\zeta}_{j}(t)$. Notice the
extra translation term compared to [8, equation (10)].

######  Remark 2.

The formula (29) upper bounds the support function of integrator reach set
resulting from the same initial condition and arbitrary compact
$\mathcal{U}\subset\mathbb{R}^{m}$ with
$\alpha_{j}:=\underset{\bm{u}\in\mathcal{U}}{\min}\;u_{j},\quad\beta_{j}:=\underset{\bm{u}\in\mathcal{U}}{{\max}\;u_{j}}$.
In other words, (20) will over-approximate the reach set associated with the
input set $\mathcal{U}$.

##  III Boundary of the Integrator Reach Set

In this Section, we build on the developments in Sec. I and II, to derive the
parametric representation of
$\bm{x}^{\text{bdy}}\in\partial\mathcal{R}\left(\\{\bm{x}_{0}\\},t\right)$,
the boundary of the forward reach set for integrator dynamics with box-valued
uncertainties in control inputs. We then provide a principled approach to
deduce the implicit equations of the boundary.

###  III-A Parameterization of the Boundary

To derive the parametric equations of the boundary of the integrator reach set
(15), we start with [3, p. 111]

$\displaystyle
h_{\mathcal{R}\left(\\{\bm{x}_{0}\\},t\right)}^{*}\left(\bm{x}^{\textup{bdy}}\right)=0,$
(35)

where $h_{\mathcal{R}\left(\\{\bm{x}_{0}\\},t\right)}^{*}$ denotes the
Legendre-Fenchel conjugate of the support function (34). That (35) yields the
equations for the boundary, follows from (19).

At any
$\bm{x}^{\textup{bdy}}\in\partial\mathcal{R}\left(\\{\bm{x}_{0}\\},t\right)$,
consider a supporting hyperplane
$\langle\bm{y},\bm{x}^{\textup{bdy}}\rangle=h_{\mathcal{R}\left(\\{\bm{x}_{0}\\},t\right)}(\bm{y})$
with outward normal $\bm{y}\in\mathbb{R}^{d}$. Using (21) and (34), the
equation (35) in our case, becomes

$\displaystyle\sum_{j=1}^{m}\>\underset{\bm{y}_{j}\in\mathbb{R}^{r_{j}}}{\inf}\\!\bigg{\\{}\\!\langle-\bm{x}^{\textup{bdy}}_{j}+\exp\left(t\bm{A}_{j}\right)\bm{x}_{j0}+\nu_{j}\bm{\zeta}_{j}(t),\bm{y}_{j}\rangle$
$\displaystyle\qquad\qquad\qquad\qquad\qquad+\mu_{j}\\!\\!\int_{0}^{t}\\!\\!|\langle\bm{y}_{j},\bm{\xi}_{j}(s)\rangle|\>{\rm{d}}s\\!\bigg{\\}}=0.$
(36)

In the Theorem 2 that follows, we use (36) to derive the parametric equations
for
$\bm{x}^{\text{bdy}}\in\partial\mathcal{R}\left(\\{\bm{x}_{0}\\},t\right)$.

######  Theorem 2.

Consider the reach set (15) with singleton
$\mathcal{X}_{0}\equiv\\{\bm{x}_{0}\\}$ and $\mathcal{U}$ given by (14). The
fixed vector $\bm{x}_{0}\in\mathbb{R}^{d}$ comprises of subvectors
$\bm{x}_{j0}\in\mathbb{R}^{r_{j}}$, where $j=1,\ldots,m$, and the relative
degree vector $\bm{r}=(r_{1},\ldots,r_{m})^{\top}\in\mathbb{N}^{m}$. Define
$\\{\mu_{j}\\}_{j=1}^{m}$ and $\\{\nu_{j}\\}_{j=1}^{m}$ as in (28). Let the
indicator function $\mathbf{1}_{k\leq\ell}:=1$ for $k\leq\ell$, and $:=0$
otherwise. Then the components of

$\bm{x}^{\textup{bdy}}=\begin{pmatrix}\bm{x}^{\textup{bdy}}_{1}\\\
\bm{x}^{\textup{bdy}}_{2}\\\ \vdots\\\
\bm{x}^{\textup{bdy}}_{m}\end{pmatrix}\in\partial\mathcal{R}\left(\\{\bm{x}_{0}\\},t\right),\;\bm{x}^{\textup{bdy}}_{j}\in\mathbb{R}^{r_{j}},\;j=1,\ldots,m,$

admit parameterization

$\displaystyle\bm{x}^{\textup{bdy}}_{j}\left(k\right)=\sum_{\ell=1}^{r_{j}}\mathbf{1}_{k\leq\ell}\>\frac{t^{\ell-k}}{(\ell-k)!}\>\bm{x}_{j0}(\ell)+\frac{\nu_{j}\>t^{r_{j}-k+1}}{(r_{j}-k+1)!}$
$\displaystyle\pm\frac{\mu_{j}}{(r_{j}-k+1)!}\bigg{\\{}(-1)^{r_{j}-1}\>t^{r_{j}-k+1}+2\sum_{q=1}^{r_{j}-1}\\!\\!(-1)^{q+1}\>s_{q}^{r_{j}-k+1}\bigg{\\}},$
(37)

where $\bm{x}^{\textup{bdy}}_{j}\left(k\right)$ denotes the $k$th component of
the $j$th subvector $\bm{x}^{\textup{bdy}}_{j}$ for $k=1,\ldots,r_{j}$. The
parameters $(s_{1},s_{2},\ldots,s_{r_{j}-1})$ satisfy $0\leq s_{1}\leq
s_{2}\leq\ldots\leq s_{r_{j}-1}\leq t$.

###### Proof.

Each of the $m$ objectives in (36) contains the absolute value of a polynomial
in $s$ that can have at most $r_{j}-1$ roots and therefore at most $r_{j}-1$
sign changes in the interval $[0,t]$. We denote these roots as $s_{1}\leq
s_{2}\leq\ldots\leq s_{r_{j}-1}$, and write:

$\displaystyle\int_{0}^{t}\\!\\!\\!\\!|\langle\bm{y}_{j},\bm{\xi}_{j}(s)\rangle|\>{\rm{d}}s=\pm\\!\int_{0}^{s_{1}}\\!\\!\\!\\!\langle\bm{y}_{j},\bm{\xi}_{j}(s)\rangle\>{\rm{d}}s\mp\\!\int_{s_{1}}^{s_{2}}\\!\\!\\!\\!\langle\bm{y}_{j},\bm{\xi}_{j}(s)\rangle\>{\rm{d}}s$
$\displaystyle\qquad\qquad\qquad\quad\quad\quad~{}~{}\pm\ldots\pm(-1)^{r_{j}-1}\\!\int_{s_{r_{j}-1}}^{t}\\!\\!\\!\\!\\!\\!\\!\\!\\!\langle\bm{y}_{j},\bm{\xi}_{j}(s)\rangle\>{\rm{d}}s$
$\displaystyle=\langle\bm{y}_{j},\pm\bm{\zeta}_{j}(0,s_{1})\mp\bm{\zeta}_{j}(s_{1},s_{2})\pm\ldots\pm(-1)^{r_{j}-1}\bm{\zeta}_{j}(s_{r_{j}-1},t)\rangle.$
(38)

Expression (38) holds even if the number of roots in $[0,t]$ is strictly less
than $r_{j}-1$ as in that case, the corresponding summand integrals become
zero.

Combining (38) and (36) we obtain

$\displaystyle\sum_{j=1}^{m}\>\underset{\bm{y}_{j}\in\mathbb{R}^{r_{j}}}{\inf}\langle-\bm{x}^{\textup{bdy}}_{j}+\exp\left(t\bm{A}_{j}\right)\bm{x}_{j0}+\nu_{j}\bm{\zeta}_{j}(t)\pm\mu_{j}\bm{\zeta}_{j}(0,s_{1})$
$\displaystyle\mp\mu_{j}\bm{\zeta}_{j}(s_{1},s_{2})\pm\ldots\pm\mu_{j}(-1)^{r_{j}-1}\bm{\zeta}_{j}(s_{r_{j}-1},t),\bm{y}_{j}\rangle=0.$
(39)

The left-hand-side of (39) is the summation of infimum values of linear
functions. Therefore, the right-hand-side of (39) can achieve zero if and only
if each of those linear function equals to zero, i.e., if and only if

$\displaystyle\bm{x}^{\textup{bdy}}_{j}=\exp\left(t\bm{A}_{j}\right)\bm{x}_{j0}+\nu_{j}\bm{\zeta}_{j}(t)\pm\mu_{j}\bm{\zeta}_{j}(0,s_{1})\mp\mu_{j}\bm{\zeta}_{j}(s_{1},s_{2})$
$\displaystyle\qquad\quad\pm\ldots\pm(-1)^{r_{j}-1}\mu_{j}\bm{\zeta}_{j}(s_{r_{j}-1},t).$
(40)

Using (2), (16) and (27), we simplify (40) to (37). This completes the proof.
∎

The following is an immediate consequence of Theorem 2.

######  Corollary 3.

The single input integrator reach set
$\mathcal{R}_{j}\left(\\{\bm{x}_{0}\\},t\right)\subset\mathbb{R}^{r_{j}}$ has
two bounding surfaces for each $j=1,\ldots,m$. In other words, there exist
$p_{j}^{\textup{upper}},p_{j}^{\textup{lower}}:\mathbb{R}^{r_{j}}\mapsto\mathbb{R}$
such that

$\mathcal{R}_{j}\left(\\{\bm{x}_{0}\\},t\right)=\\{\bm{x}\in\mathbb{R}^{r_{j}}\mid
p_{j}^{\textup{upper}}(\bm{x})\leq 0,\;p_{j}^{\textup{lower}}(\bm{x})\leq
0\\},$

with boundary
$\partial\mathcal{R}_{j}\left(\\{\bm{x}_{0}\\},t\right)=\\{\bm{x}\in\mathbb{R}^{r_{j}}\mid
p_{j}^{\textup{upper}}(\bm{x})=0\\}\cup\\{\bm{x}\in\mathbb{R}^{r_{j}}\mid
p_{j}^{\textup{lower}}(\bm{x})=0\\}$.

###### Proof.

From (37), we get two different parametric representations of
$\bm{x}_{j}^{\text{bdy}}$ in terms of $(s_{1},s_{2},\ldots,s_{r_{j}-1})$. One
parametric representation results from the choice of positive sign for the
$\pm$ appearing in (37), and another for the choice of negative sign for the
same. The plus (resp. minus) sign recovers $p_{j}^{\textup{upper}}(\bm{x})=0$
(resp. $p_{j}^{\textup{lower}}(\bm{x})=0$) in parametric representation (37).
∎

Let us illustrate the boundary parametrization (37) for the case
$\bm{r}=(r_{1},r_{2})^{\top}=(2,3)^{\top}$. In this case, we get

Figure 1: The integrator reach set
$\mathcal{R}(\\{\bm{x}_{0}\\},t)\subset\mathbb{R}^{2}$ with $d=2$, $m=1$,
$\bm{x}_{0}=(0.05,0.05)^{\top}$, $\mathcal{U}\equiv[\alpha,\beta]=[-1,1]$ at
$t=2.1$. Figure 2: The integrator reach set
$\mathcal{R}(\\{\bm{x}_{0}\\},t)\subset\mathbb{R}^{3}$ with $d=3$, $m=1$,
$\bm{x}_{0}=(0.1,0.2,0.3)^{\top}$, $\mathcal{U}\equiv[\alpha,\beta]=[-1,1]$ at
$t=2.1$. The wireframes correspond to the upper and lower surfaces.

$\displaystyle\\!\\!\begin{pmatrix}\bm{x}^{\text{bdy}}_{1}(1)\\\
\vspace*{-0.05in}\\\
\bm{x}^{\text{bdy}}_{1}(2)\end{pmatrix}\\!\\!=\\!\\!\begin{pmatrix}\bm{x}_{10}(1)+t\bm{x}_{10}(2)+\nu_{1}(t^{2}/2)\pm\mu_{1}\left(s_{1}^{2}-t^{2}/2\right)\\\
\vspace*{-0.05in}\\\
\bm{x}_{10}(2)+\nu_{1}t\pm\mu_{1}\left(2s_{1}-t\right)\end{pmatrix},$ (41)

and

$\displaystyle\\!\\!\begin{pmatrix}\bm{x}^{\text{bdy}}_{2}(1)\\\
\vspace*{-0.05in}\\\ \bm{x}^{\text{bdy}}_{2}(2)\\\ \vspace*{-0.05in}\\\
\bm{x}^{\text{bdy}}_{2}(3)\end{pmatrix}=\begin{pmatrix}\bm{x}_{20}(1)+t\bm{x}_{20}(2)+(t^{2}/2)\bm{x}_{20}(3)\\\
+\nu_{2}(t^{3}/6)\pm\mu_{2}\left(t^{3}/6+2s_{1}^{3}/6-2s_{2}^{3}/6\right)\\\
\vspace*{-0.05in}\\\
\bm{x}_{20}(2)+t\bm{x}_{20}(3)+\nu_{2}(t^{2}/2)\pm\mu_{2}\left(t^{2}/2\right.\\\
\left.+2s_{1}^{2}/2-2s_{2}^{2}/2\right)\\\ \vspace*{-0.05in}\\\
\bm{x}_{20}(3)+\nu_{2}t\pm\mu_{2}\left(t+2s_{1}-2s_{2}\right)\end{pmatrix}.$
(42)

In (41), taking plus (resp. minus) signs in each of component gives the
parametric representation of the curve $p_{1}^{\text{upper}}=0$ (resp.
$p_{1}^{\text{lower}}=0$). The union of these curves defines
$\partial\mathcal{R}_{1}$ and are depicted in Fig. 1. We note that the
parameterization (41) appeared in [3, p. 111].

Likewise, in (42), taking plus (resp. minus) signs in each of component gives
the parametric representation of the surface $p_{2}^{\text{upper}}(\bm{x})=0$
(resp. $p_{2}^{\text{lower}}=0$). The resulting set $\mathcal{R}_{2}$ is the
triple integrator reach set, and is shown in Fig. 2. The single input
integrator reach set for $\bm{r}=(r_{1},r_{2})^{\top}=(2,1)^{\top}$ is shown
in Fig. 3.

Next, we deduce the implicit equations of the boundary associated with the
parameterization (37).

###  III-B Implicitization of the Boundary

We recall the concepts of _variety_ and _ideal_ [9, Ch. 1]. Let
$p_{1},\ldots,p_{n}\in\mathbb{R}[x_{1},\ldots,x_{d}]$, the vector space of
real-valued $d$-variate polynomials. The (affine) variety
$V_{\mathbb{R}[x_{1},\ldots,x_{d}]}(p_{1},\ldots,p_{n})$ is the set of all
solutions of the system
$p_{1}(x_{1},x_{2},\ldots,x_{d})=\ldots=p_{n}(x_{1},x_{2},\ldots,x_{d})=0$.

Given $p_{1},\ldots,p_{n}\in\mathbb{R}[x_{1},\ldots,x_{d}]$, the set

$I:=\bigg{\\{}\sum_{i=1}^{n}\alpha_{i}p_{i}\mid\alpha_{1},\ldots,\alpha_{n}\in\mathbb{R}[x_{1},\ldots,x_{d}]\bigg{\\}}$

is called the ideal generated by $p_{1},\ldots,p_{n}$. We write this
symbolically as $I=\langle\langle p_{1},\ldots,p_{n}\rangle\rangle$.
Intuitively, $\langle\langle p_{1},\ldots,p_{n}\rangle\rangle$ can be thought
of as the set of all polynomial consequences of the given system of $n$
polynomial equations in $d$ indeterminates.

We notice that (37) gives polynomial parameterizations of the components of
$\bm{x}_{j}^{\text{bdy}}$ for all $j=1,\ldots,m$. In particular, for each
$k\in\\{1,\ldots,r_{j}\\}$, the-right-hand-side of (37) is a homogeneous
polynomial in $r_{j}-1$ parameters $(s_{1},s_{2},\ldots,s_{r_{j}-1})$ of
degree $r_{j}-k+1$. By polynomial implicitization [9, p. 134], the
corresponding implicit equations
$p_{j}^{\text{upper}}\left(\bm{x}_{j}^{\text{bdy}}\right)=0$ (when fixing plus
sign for $\pm$ in (37)) and
$p_{j}^{\text{lower}}\left(\bm{x}_{j}^{\text{bdy}}\right)=0$ (when fixing
minus sign for $\pm$ in (37)), must define affine varieties
$V_{\mathbb{R}[x_{1},...,x_{r_{j}}]}(p_{j}^{\text{upper}}),V_{\mathbb{R}[x_{1},...,x_{r_{j}}]}(p_{j}^{\text{lower}})$
in $\mathbb{R}\left[x_{1},\ldots,x_{d}\right]$. Theorem 4 below provides the
implicit equations associated with the parameterization (37).

Figure 3: The integrator reach set $\mathcal{R}(\\{\bm{x}_{0}\\},t)$ with
$d=3$, $m=2$, $\bm{r}=(2,1)^{\top}$, $\bm{x}_{0}=(1,1,0)^{\top}$,
$[\alpha_{1},\beta_{1}]=[-5,5]$, $[\alpha_{2},\beta_{2}]=[-3,3]$ at $t=4$.

######  Theorem 4.

Consider $d$ dimensional integrator reach set (15) resulting from singleton
$\mathcal{X}_{0}\equiv\\{\bm{x}_{0}\\}$ and $\mathcal{U}$ given by (14).
Denote the corresponding single input integrator reach sets
$\mathcal{R}_{1},\ldots,\mathcal{R}_{m}$ as in (20). Let
$\delta_{j}:=\lfloor\frac{r_{j}-1}{2}\rfloor$ and

$\displaystyle\rho_{j,k}^{\pm}:=$
$\displaystyle\frac{(r_{j}-k+1)!}{2\mu_{j}}\bigg{\\{}\bm{x}^{\textup{bdy}}_{j}\left(k\right)-\sum_{\ell=1}^{r_{j}}\mathbf{1}_{k\leq\ell}\>\frac{t^{\ell-k}}{(\ell-k)!}\>\bm{x}_{j0}(\ell)\bigg{\\}}$
$\displaystyle-\frac{1}{2}\bigg{\\{}\pm(-1)^{r_{j}-1}\>t^{r_{j}-k+1}+\frac{\nu_{j}}{\mu_{j}}\>t^{r_{j}-k+1}\bigg{\\}},$
(43)

where $\bm{r}=(r_{1},\ldots,r_{m})^{\top}\in\mathbb{N}^{m}$ is the relative
degree vector, and $\bm{x}_{0}\in\mathbb{R}^{d}$ comprises of subvectors
$\bm{x}_{j0}\in\mathbb{R}^{r_{j}}$, $j=1,\ldots,m$. The algebraic
hypersurfaces
$p_{j}^{\text{upper}},p_{j}^{\text{lower}}\in\mathbb{R}[x_{1},...,x_{r_{j}}]$
bounding $\mathcal{R}_{j}\left(\\{\bm{x}_{0}\\},t\right)$ have degree
$\left(\delta_{j}+1\right)\left(r_{j}-\delta_{j}\right)$, and are given by the
following vanishing Hankel determinant

$\displaystyle{\mathrm{det}}[A_{j}^{\pm}(r_{j}-2\delta_{j}+i+j)]_{i,j=0}^{\delta_{j}}=0,$
(44)

where $A^{\pm}_{j}$’s are obtained from

$\displaystyle\sum_{k\geq
0}A^{\pm}_{j}(k)\tau^{k}=\exp\left(-\sum_{k=1}^{r_{j}}\frac{\lambda_{j}^{\pm}(k)}{k}\tau^{k}\right),$
(45)

and $\lambda_{j}^{\pm}(k):=\rho_{r_{j}-k+1}^{\pm}$ for $k=1,\ldots,r_{j}$.

###### Proof.

To ease notation, let us drop the superscript $\pm$ from
$\rho_{r_{j}-k+1}^{\pm},A^{\pm}_{j}$ and $\lambda_{j}^{\pm}(k)$ in (43). These
superscripts can be added back at the end of the proof, without loss of
generality.

Following [10], we next define the sequence $A_{j}(k)$ (the sequence is with
respect to index $k$) in terms of $s_{1},s_{2},\dots,s_{r_{j}-1}$ via the
generating function (see e.g.,[11, Ch. 1])

$\displaystyle F_{j}(\tau)=\sum_{k\geq
0}A_{j}(k)\tau^{k}=\frac{(1-s_{1}\tau)(1-s_{3}\tau)\cdots}{(1-s_{2}\tau)(1-s_{4}\tau)\cdots}.$
(46)

Taking the logarithmic derivative of (46), and using the generating functions
$(1-s_{q}\tau)^{-1}=\sum_{k\geq 0}\left(s_{q}\tau\right)^{k}$ for all
$q=1,\ldots,r_{j}-1$, yields

$\displaystyle\\!\\!\frac{F_{j}^{\prime}(\tau)}{F_{j}(\tau)}=-s_{1}\\!\sum_{k\geq
0}\\!\left(s_{1}\tau\right)^{k}\\!+s_{2}\\!\sum_{k\geq
0}\\!\left(s_{2}\tau\right)^{k}\\!-s_{3}\\!\sum_{k\geq
0}\\!\left(s_{3}\tau\right)^{k}\\!+\ldots.$ (47)

In (47), the superscript ′ denotes derivative with respect to $\tau$.
Integrating (47) with respect to $\tau$, we get

$\displaystyle
F_{j}(\tau)=\exp\left(-\sum_{k=1}^{r_{j}}\frac{\lambda_{j}(k)}{k}\tau^{k}\right).$
(48)

Equating (46) and (48), we obtain (45).

On the other hand, since the generating function (46) is a rational function
with denominator polynomial of degree $\delta_{j}$, the Hankel determinant in
(44) vanishes[12],[13, p. 5, Lemma III]. Finally, reverting back the
$\lambda_{j}$’s to the $\rho_{j}$’s, and inserting back the plus-minus
superscripts, result in the desired implicit polynomial, say
$\wp(\rho_{j,1}^{\pm},\rho_{j,2}^{\pm},\ldots,\rho_{j,r_{j}}^{\pm})$, of
degree $(\delta_{j}+1)(r_{j}-\delta_{j})$. ∎

To illustrate Theorem 4, let $d=3$ and $m=1$; the latter allows us to drop the
subscript $j$. Then (44) becomes

$\displaystyle{\mathrm{det}}\left(\begin{bmatrix}A^{\pm}(1)&A^{\pm}(2)\\\
A^{\pm}(2)&A^{\pm}(3)\end{bmatrix}\right)=0.$ (49)

In this case, equating (46) and (48) gives

$\displaystyle
A^{\pm}(1)=-\lambda^{\pm}(1),\;A^{\pm}(2)=\frac{1}{2}\lambda^{\pm}({1})^{2}-\frac{1}{2}\lambda^{\pm}({2}),\;$
$\displaystyle
A^{\pm}({3})=-\frac{1}{6}\lambda^{\pm}({1})^{3}+\frac{1}{2}\lambda^{\pm}({1})\lambda^{\pm}({2})-\frac{1}{3}\lambda^{\pm}({3}).$

Substituting the above back in (49) yields the quartic polynomial
$\lambda^{\pm}({1})^{4}+3\lambda^{\pm}({2})^{2}-4\lambda^{\pm}({3})\lambda^{\pm}({1})=0$,
which under the mapping
$(\lambda^{\pm}(1),\lambda^{\pm}(2),\lambda^{\pm}(3))\mapsto(\rho^{\pm}(3),\rho^{\pm}(2),\rho^{\pm}(1))$,
yields degree 4 implicitized polynomial

$\displaystyle\wp(\rho^{\pm}(1),\rho^{\pm}(2),\rho^{\pm}(3))$
$\displaystyle=\rho^{\pm}({3})^{4}-4\rho^{\pm}({3})\rho^{\pm}({1})+3\rho^{\pm}({2})^{2}=0.$
(50)

For $k=1,2,3$, substituting for the
$\rho^{\pm}({1}),\rho^{\pm}({2}),\rho^{\pm}({3})$ in (50) from (43), allows us
to recover $p_{j}^{\text{upper}}$ and $p_{j}^{\text{lower}}$.

We can verify (50) using elementary algebra by eliminating the parameters
$\left(s_{1},s_{2}\right)$ from (42). However, it is practically impossible to
derive the implicitization via brute force algebra for $d=4$ or higher.

In summary, (44) gives the implicitization of the bounding hypersurfaces of
the single input integrator reach set (up to the change of variables). The
implicitization of the hypersurfaces in multi-input case, then, is the
Cartesian product of these single input implicit hypersurfaces.

##  IV Taxonomy of the Integrator Reach Set

In this Section, we address the question: what type of compact convex set the
integrator reach set resulting from (14) is? This question is particularly
appealing since several types of convex sets are well-studied in the convex
analysis literature. In systems-control literature too, the knowledge of the
type of convex set (e.g., semialgebraic, spectrahedron) under investigation,
is often leveraged to facilitate analysis, and to design specialized
algorithms.

We recall the concept of semialgebraic sets. A set in $\mathbb{R}^{d}$ is
called _basic semialgebraic_ if it can be written as a finite conjunction of
polynomial inequalities and equalities, the polynomials being in
$\mathbb{R}\left[x_{1},\ldots,x_{d}\right]$. Finite union of basic
semialgebraic sets is called a _semialgebraic set_. A semialgebraic set need
not be basic semialgebriac; see e.g., [14, Example 2.2]. Our first taxonomy
result is the following.

######  Theorem 5.

The integrator reach set (15) with $\mathcal{X}_{0}\equiv\\{\bm{x}_{0}\\}$ and
box-valued input set $\mathcal{U}$, is semialgebraic at any time $t$.

###### Proof.

The implicitization derived in Sec. III-B, by construction, demonstrated that
$p_{j}^{\text{upper}},p_{j}^{\text{lower}}\in\mathbb{R}\left[x_{1},\ldots,x_{r_{j}}\right]$
for all $j=1,\ldots,m$. In words, $p_{j}^{\text{upper}},p_{j}^{\text{lower}}$
are real algebraic hypersurfaces for all $j=1,\ldots,m$. For a different proof
that uses the parametric, instead of the implicit representation of the
boundary, see [15, Appendix E]. Alternatively, the statement can be proved by
applying Tarski-Seidenberg theorem [16, Ch. 1] on (37). ∎

So far, we have established that the compact convex set
$\mathcal{R}\left(\\{\bm{x}_{0}\\},t\right)$ is semialgebraic (Theorem 5), and
a translated zonoid (Remark 1). Convex semialgebraic sets contain a well-known
subclass called _spectrahedra_ , also referred to as _linear matrix inequality
(LMI) representable sets_. Spectrahedra are affine slices of the symmetric
positive semidefinite cone.

The projections of spectrahedra are referred to as _spectrahedral shadows_ ,
which subsume the class of spectrahedra. It is possible that a compact convex
set is spectrahedral shadow, i.e., admits LMI representation when lifted to a
higher dimensional space, but on its own, not a spectrahedron in the original
dimensions. For an example of spectrahedral shadow that is not a
spectrahedron, see [17, p. 23].

Helton and Vinnikov [18, Thm. 3.1] showed that for a $d$ dimensional compact
convex set to be spectrahedron, a necessary condition is _rigid convexity_. A
compact semialgebriac set is rigidly convex if the number of intersections
made by a generic line passing through an interior point of the set, with its
real algebraic boundary is equal to the degree of the bounding algebraic
hypersurfaces; see [18, Sec. 3.1 and 3.2].

We know from Sec. III-B that the $d$ dimensional single input (i.e., $m=1$)
integrator reach set has bounding algebraic surfaces of degree
$(\lfloor\frac{d-1}{2}\rfloor+1)(d-\lfloor\frac{d-1}{2}\rfloor)$. In
particular, for $d=2$, Fig. 4(a) shows that a generic line has $4$
intersections with the bounding real algebraic curves whereas from Theorem 4,
we know that $p^{\text{upper}},p^{\text{lower}}$ in this case, are degree $2$
polynomials.

Likewise, for $d=3$ and $m=1$, Fig. 4(b) reveals that a generic line has $6$
intersections with the bounding real algebraic surfaces whereas from (50), we
know that the polynomials $p^{\text{upper}},p^{\text{lower}}$ in this case,
are of degree $4$. So the integrator reach set box-valued input set
$\mathcal{U}$, is not a spectrahedron.

(a) Real algebraic curves $p^{\text{upper}},\allowbreak p^{\text{lower}}$ for
the double integrator.

(b) Real algebraic surfaces $p^{\text{upper}},\allowbreak p^{\text{lower}}$
for the triple integrator.

Figure 4: The bounding polynomials for the double and triple integrator reach
sets at $t=0.5$ with $\bm{x}_{0}=\bm{0}$ and $\mu=1$. Figure 5: The summary of
the taxonomy for the integrator reach set with box-valued input set
$\mathcal{U}$.

An interesting question is the following: could this reach set be
spectrahedral shadow? Although some calculations show that _sufficient_
conditions as in [17] do not seem to hold, the authors do not have conclusive
answer to this matter. The difficulty here, to the best of the authors’
knowledge, is that the necessary conditions for a compact semialgebraic set to
be spectrahedral shadow, is not available in the current literature. This
makes it hard to falsify this possibility for the integrator reach set.

We graphically summarize our taxonomy results in Fig. 5; the highlighted
region shows where the integrator reach set belongs. Our future work will
explore the possibility of further narrowing down the highlighted region.

##  V Conclusions and Future Work

The present paper is part of a research program [8, 15] to understand the
specific geometry of the integrator reach sets. This is motivated by the fact
that integrator reach sets feature prominently for benchmarking the
performance of reach set over-approximation algorithms, and that they may be
used to over-approximate the forward reach sets of differentially flat
nonlinear systems in the normal coordinates (and then numerically map the
boundary back in original coorodinates via known coordinate transforms).

In the present paper, we deduced a closed-form formula for the support
function of the forward reach set of the integrator dynamics with box-valued
uncertainties in the control inputs. Using the support function formula, we
then derived analytic expressions of the boundary of these compact convex sets
in both parametric and implicit forms. We then delved into the taxonomy
question. We argued that this reach set is a translated zonoid, semialgebraic,
and not a spectrahedron.

Our future work will investigate the possibility of further narrowing down the
highlighted region in Fig. 5. In addition, using the ideas presented herein to
design algorithms for computing or tightly over-approximating the forward
reach sets of differentially flat systems, is a topic of the authors’ ongoing
research.

## References

  * [1] D. Maksarov and J. Norton, “State bounding with ellipsoidal set description of the uncertainty,” _International Journal of Control_ , vol. 65, no. 5, pp. 847–866, 1996.
  * [2] N. Kochdumpe, V. Gaßmann, F. Gruber, M. Klischat, B. Schürman, and M. Althoff, “AROC 2020 manual,” Technical University of Munich. [Online]. Available: https://tumcps.github.io/AROC/data/Aroc2020Manual.pdf
  * [3] A. B. Kurzhanski and P. Varaiya, _Dynamics and Control of Trajectory Tubes: Theory and Computation_. Springer, 2014, vol. 85.
  * [4] J. Lévine, P. Martin, and P. Rouchon, “Flatness and defect of non-linear systems: introductory theory and examples,” _International Journal of Control_ , vol. 61, no. 6, pp. 1327–1361, 1995.
  * [5] R. M. Murray, M. Rathinam, and W. Sluis, “Differential flatness of mechanical control systems: a catalog of prototype systems.” ASME International Mechanical Engineering Congress and Exposition, 1995.
  * [6] J. Lévine, P. Martin, and P. Rouchon, “A Lie-Bäcklund approach to equivalence and flatness of nonlinear systems,” _IEEE Transactions on Automatic Control_ , vol. 44, no. 5, pp. 922–937, 1999.
  * [7] P. Varaiya, “Reach set computation using optimal control,” in _Verification of Digital and Hybrid Systems_. Springer, 2000, pp. 323–331.
  * [8] S. Haddad and A. Halder, “The convex geometry of integrator reach sets,” in _2020 American Control Conference (ACC)_. IEEE, 2020, pp. 4466–4471.
  * [9] D. Cox, J. Little, and D. OShea, _Ideals, varieties, and algorithms: an introduction to computational algebraic geometry and commutative algebra_. Springer Science & Business Media, 2013.
  * [10] G. Zaimi, “A polynomial implicitization,” MathOverflow. [Online]. Available: https://mathoverflow.net/q/381335
  * [11] H. S. Wilf, _generatingfunctionology_. CRC press, 2005.
  * [12] L. Kronecker, _Zur Theorie der Elimination einer Variabeln aus zwei algebraischen Gleichungen_. Buchdruckerei der Königl. Akademie der Wissenschaften (G. Vogt), 1881.
  * [13] R. Salem, _Algebraic numbers and Fourier analysis_. Wadsworth Publishing Company, 1983.
  * [14] C. Vinzant, “The geometry of spectrahedra,” in _Sum of Squares: Theory and Applications_. American Mathematical Society, 2020, pp. 11–36.
  * [15] S. Haddad and A. Halder, “The curious case of integrator reach sets, Part I: Basic theory,” _arXiv preprint, arXiv:eess.SY/2102.11423_ , 2021. [Online]. Available: https://arxiv.org/pdf/2102.11423.pdf
  * [16] J. Bochnak, M. Coste, and M.-F. Roy, _Real algebraic geometry_. Springer Science & Business Media, 2013, vol. 36.
  * [17] J. W. Helton and J. Nie, “Semidefinite representation of convex sets,” _Mathematical Programming_ , vol. 122, no. 1, pp. 21–64, 2010.
  * [18] J. W. Helton and V. Vinnikov, “Linear matrix inequality representation of sets,” _Communications on Pure and Applied Mathematics_ , vol. 60, no. 5, pp. 654–674, 2007.