prediction analogies can only be drawn in relation to ideal-gas transient-
growth results.
In the non-transcritical cases, transient growth is marginally affected by the
non-ideal thermodynamic regime. In the best-case-scenarios (subcritical wall-
heating and supercritical wall-cooling), the levels of amplification at which
modal growth prevails over transient growth are around $N\approx 4$.
In the transcritical regime, due to the highly inflectional base flow, modal
instability dominates, especially for transcritical wall-cooling, where no
subcritical transition is likely to occur. In this case, the increase in the
non-modal $N$-factor, i.e., $N_{opt}$, is significant compared to the non-
transcritical wall-cooling cases, sharing the same effect as an adverse-
pressure gradient in ideal-gas transient-growth calculations. A similar
analogy can be drawn for transcritical wall-heating, but with lower transient-
growth amplifications, reaching up to $N\approx 3.5$ before a strong modal
amplification emerges. Nevertheless, a critical Reynolds number of
$Re_{x}\approx 3\times 10^{5}$ could be sufficiently large to make transition
via transient growth the most likely transition scenario when the level of
free-stream turbulence is high.
The increase in supercritical pressure does not affect the location of the
generalized inflection point. However, a rise in $p$ causes the onset of modal
instability to shift towards higher Reynolds numbers. For instance, in the
transcritical wall-heating, at
$p^{*}=$85\text{\,}\mathrm{b}\mathrm{a}\mathrm{r}$$, Mode I is the most
unstable mode, contrary to the usual Mode II. Once more, $N_{opt.}$ is hardly
influenced by a property change, specifically pressure. The Reynolds number at
which non-modal growth is overtaken by modal growth increases by a factor of
$3$ from $p^{*}=$80\text{\,}\mathrm{b}\mathrm{a}\mathrm{r}$$ to
$85\text{\,}\mathrm{b}\mathrm{a}\mathrm{r}$, enhancing the chances of
subcritical transition for transcritical wall-heating due to strong initial
disturbances. Conversely, in the case of transcritical wall-cooling, modal
growth always dominates.

###### Acknowledgements.

The European Research Council (grant no. ERC-2019-CoG-864660, Critical) has
supported the investigations presented in this paper. The authors acknowledge
the use of computational resources of DelftBlue supercomputer, provided by
Delft High Performance Computing Centre (https://www.tudelft.nl/dhpc).

##  Appendix A Supercritical $\text{CO}_{2}$

The material dependent parameters of $\text{CO}_{2}$ are provided in Table 5.

Table 5: Thermodynamic properties of $\text{CO}_{2}$: gas constant ($R^{*}_{g}$), critical pressure ($p^{*}_{c}$), and critical temperature ($T^{*}_{c}$). $R^{*}_{g}\ ($\text{\,}\mathrm{J}\mathrm{/}\mathrm{k}\mathrm{g}\mathrm{/}\mathrm{K}$)$ | $p^{*}_{c}\ ($\text{\,}\mathrm{b}\mathrm{a}\mathrm{r}$)$ | $T^{*}_{c}\ ($\text{\,}\mathrm{K}$)$
---|---|---
$188.9$ | $73.9$ | $304.1$

The evolution of the pseudo-boiling temperature as a function of pressure and
temperature is plotted in Fig. 29 along with contours of $c^{*}_{p}$.

Figure 29: Evolution of the pseudo-boiling temperature $T^{*}_{pc}$ ( ) over
$p^{*}$ with contours of $c^{*}_{p}$. Three isobars ( ) are drawn: (a)
$T^{*}_{pc}(p^{*}_{c}=$80\text{\,}\mathrm{b}\mathrm{a}\mathrm{r}$)=$307.7\text{\,}\mathrm{K}$$,
(b)
$T^{*}_{pc}(p^{*}_{c}=$81.15\text{\,}\mathrm{b}\mathrm{a}\mathrm{r}$)=$308.3\text{\,}\mathrm{K}$$,
and (c)
$T^{*}_{pc}(p^{*}_{c}=$84.8\text{\,}\mathrm{b}\mathrm{a}\mathrm{r}$)=$310.4\text{\,}\mathrm{K}$$.

##  Appendix B Modal stability: base-flow matrices

The adjusted terms in the base-flow matrices of Ref. [7] for the stability
equation (see Eq. (10)) are reported hereafter. For simplicity, the derivative
of a thermodynamic quantity with respect to $\bar{p}$ at constant $\bar{T}$,
and vice versa, is denoted as $\partial/\partial\bar{p}$ instead of
$\partial/\partial\bar{p}|_{\bar{T}}$ and as $\partial/\partial\bar{T}$
instead of $\partial/\partial\bar{T}|_{\bar{p}}$, respectively. Furthermore,
any occurrence of $\partial/\partial\bar{\rho}$ in the matrices of Ref. [7]
changes to $\partial/\partial\bar{p}$ and will not be reported hereafter. The
modified elements are:

$\left.\begin{array}[]{cc}\mathbf{L_{t}}(1,1)=\dfrac{\partial\bar{\rho}}{\partial\bar{p}},\quad\mathbf{L_{t}}(1,5)=\dfrac{\partial\bar{T}}{\partial\bar{p}},\end{array}\right\\}$
(33)

$\left.\begin{array}[]{cc}\mathbf{L_{x}}(1,1)=\bar{u}\dfrac{\partial\bar{\rho}}{\partial\bar{p}},\quad\mathbf{L_{x}}(1,5)=\bar{u}\dfrac{\partial\bar{T}}{\partial\bar{p}},\\\\[8.61108pt]
\mathbf{L_{x}}(2,1)=1,\quad\mathbf{L_{x}}(2,5)=0,\\\\[8.61108pt]
\mathbf{L_{x}}(5,2)=Ec_{\infty}\bar{p},\quad\mathbf{L_{x}}(5,3)=-\dfrac{2Ec_{\infty}\bar{\mu}}{Re_{\delta}}\dfrac{\partial\bar{u}}{\partial
y},\\\ \end{array}\right\\}$ (34)

$\left.\begin{array}[]{cc}\mathbf{L_{y}}(3,1)=1,\quad\mathbf{L_{y}}(3,5)=0,\\\\[8.61108pt]
\mathbf{L_{y}}(5,1)=-\dfrac{1}{Re_{\delta}Pr_{\infty}}\dfrac{\partial\bar{\kappa}}{\partial\bar{p}}\dfrac{\partial\bar{T}}{\partial
y},\quad\mathbf{L_{y}}(5,2)=-\dfrac{2Ec_{\infty}\bar{\mu}}{Re_{\delta}}\dfrac{\partial\bar{u}}{\partial
y},\\\\[8.61108pt]
\quad\mathbf{L_{y}}(5,3)=Ec_{\infty}\bar{p},\quad\mathbf{L_{y}}(5,5)=-\dfrac{1}{Re_{\delta}Pr_{\infty}}\left(\dfrac{\partial\bar{\kappa}}{\partial
y}+\dfrac{\partial\bar{\kappa}}{\partial\bar{T}}\dfrac{\partial\bar{T}}{\partial
y}\right),\\\ \end{array}\right\\}$ (35)

$\left.\begin{array}[]{cc}\mathbf{L_{z}}(4,1)=1,\quad\mathbf{L_{z}}(4,5)=0,\\\\[8.61108pt]
\mathbf{L_{z}}(5,4)=Ec_{\infty}\bar{p},\\\ \end{array}\right\\}$ (36)

$\left.\begin{array}[]{cc}\mathbf{L_{q^{\prime}}}(2,1)=-\dfrac{1}{Re_{\delta}}\left(\dfrac{\partial\bar{\mu}}{\partial\bar{p}}\dfrac{\partial^{2}\bar{u}}{\partial
y^{2}}+\dfrac{\partial^{2}\bar{\mu}}{\partial\bar{p}\,\partial\bar{T}}\dfrac{\partial\bar{T}}{\partial
y}\dfrac{\partial\bar{u}}{\partial
y}\right),\quad\mathbf{L_{q^{\prime}}}(2,5)=-\dfrac{1}{Re_{\delta}}\left(\dfrac{\partial\bar{\mu}}{\partial\bar{T}}\dfrac{\partial^{2}\bar{u}}{\partial
y^{2}}+\dfrac{\partial^{2}\bar{\mu}}{\partial\bar{T}^{2}}\dfrac{\partial\bar{T}}{\partial
y}\dfrac{\partial\bar{u}}{\partial y}\right),\\\\[12.91663pt]
\mathbf{L_{q^{\prime}}}(3,1)=0,\quad\mathbf{L_{q^{\prime}}}(3,5)=0,\\\\[8.61108pt]
\mathbf{L_{q^{\prime}}}(5,1)=-\dfrac{Ec_{\infty}}{Re_{\delta}}\left[\dfrac{\partial\bar{\mu}}{\partial\bar{p}}\left(\dfrac{\partial\bar{u}}{\partial
y}\right)^{2}+\dfrac{1}{Ec_{\infty}Pr_{\infty}}\dfrac{\partial^{2}\bar{\kappa}}{\partial\bar{p}\,\partial\bar{T}}\left(\dfrac{\partial\bar{T}}{\partial
y}\right)^{2}+\dfrac{1}{Ec_{\infty}Pr_{\infty}}\dfrac{\partial\bar{\kappa}}{\partial\bar{p}}\dfrac{\partial^{2}\bar{T}}{\partial
y^{2}}\right],\\\\[12.91663pt]
\mathbf{L_{q^{\prime}}}(5,5)=-\dfrac{Ec_{\infty}}{Re_{\delta}}\left[\dfrac{\partial\bar{\mu}}{\partial\bar{T}}\left(\dfrac{\partial\bar{u}}{\partial
y}\right)^{2}+\dfrac{1}{Ec_{\infty}Pr_{\infty}}\dfrac{\partial^{2}\bar{\kappa}}{\partial\bar{T}^{2}}\left(\dfrac{\partial\bar{T}}{\partial
y}\right)^{2}+\dfrac{1}{Ec_{\infty}Pr_{\infty}}\dfrac{\partial\bar{\kappa}}{\partial\bar{T}}\dfrac{\partial^{2}\bar{T}}{\partial
y^{2}}\right],\\\\[8.61108pt] \end{array}\right\\}$ (37)

$\left.\begin{array}[]{cc}\mathbf{V_{xx}}(5,5)=\mathbf{V_{yy}}(5,5)=\mathbf{V_{zz}}(5,5)=-\dfrac{\bar{\kappa}}{Re_{\delta}Pr_{\infty}}.\\\\[8.61108pt]
\end{array}\right\\}$ (38)

##  Appendix C Validation: transient growth

The transient-growth analysis presented in this work is validated by
reproducing the spatial results of Ref. [22] for an ideal gas (adiabatic wall,
$T^{*}_{tot}=$333\text{\,}\mathrm{K}$$, $Pr_{\infty}=0.70$,
$c^{*}_{p}/c^{*}_{v}=1.40$, $\omega=0$, $Re_{\delta}=300$). A good agreement
is observed in Fig. 30(a) at different Mach numbers.
To obtain a grid-independent solution when using a non-ideal gas, the
influence of mesh size $N$ and mesh height $y_{max}$ is investigated. A
transcritical wall-heating reference case at constant supercritical pressure
is chosen: $p^{*}_{\infty}=$80\text{\,}\mathrm{b}\mathrm{a}\mathrm{r}$$,
$T^{*}_{\infty}/T^{*}_{pc}=0.90$, isothermal wall with
$T^{*}_{w}/T^{*}_{pc}=1.05$, $M_{\infty}=10^{-3}$, $Pr_{\infty}=2.11$,
$\omega=0.013$, $\beta=0.45$, and $Re_{\delta}=300$. The results, shown in
Fig. 30(b) for $N$ and Fig. 30(c) for $y_{max}$, respectively, reveal that a
value of $N=300$ and $y_{max}=75$ is necessary for mesh independence. Hence,
this numerical setup is retained throughout the non-modal calculations in this
study.

Figure 30: (a) ideal gas, comparison of maximum energy amplification:
$M_{\infty}=0.5$ ( ), $M_{\infty}=1.0$ ( ), $M_{\infty}=2.0$ ( ),
$M_{\infty}=3.0$ ( ), Ref. [22] marked by a colored circle ($\circ$); (b) non-
ideal gas, influence of mesh size $N$ on the energy amplification $G(x)$:
$N=200$ ( ), $N=300$ ( ), $N=400$ ( ); (c) non-ideal gas, influence of mesh
height $y_{max}$ on the energy amplification $G(x)$: $y_{max}=75$ ( ),
$y_{max}=100$ ( ), $y_{max}=150$ ( ).

##  Appendix D Temporal vs. spatial transient growth

So far, spatial optimal growth has been analyzed for its simpler
interpretation in experimental facilities [22]. However, it is noteworthy to
investigate whether similar non-modal behavior is obtained regardless of the
considered optimization procedure. For the temporal analysis, initial
disturbances are optimized at $t=0$ in order to achieve the maximum temporal
energy growth at a later time $t_{max}$ over $\beta$ and streamwise wavenumber
$\alpha$ [17]. The case selected for the comparison is the transcritical case
T09w105, where oblique disturbances are optimal.
Figs. 31(a,b) and Figs. 31(c,d) display the maximum energy amplification
$G_{max}$ and the maximum location $x_{max}$ and time $t_{max}$, respectively,
for both the spatial and temporal analysis of case T09w105.

Figure 31: Transcritical case T09w105: (a) contour plots of
$G_{max}(\omega,\beta)$, (b) contour plots of $G_{max}(\alpha,\beta)$, (c)
contour plots of $x_{max}(\omega,\beta)$, and (d) contour plots of
$t_{max}(\alpha,\beta)$. $G_{opt}$ is denoted with a black star
($\filledstar$) symbol.

The qualitative behavior of $G_{max}$ is nearly identical for low values of
$\omega$ and $\alpha$ for both Fig. 31(a) and Fig. 31(b). At large $\omega$
and $\alpha$, the spatial framework reveals a larger sub-optimal growth near
the $\omega$-axis, which is negligible in the temporal framework. In the
latter, the location of $G_{opt}$, defined as
$\max\\{G_{max}(\alpha,\beta)\\}$, is found at $\beta_{opt}=0.48$ and
$\alpha_{opt}=0.014$. This confirms the observation of Ref. [23], in which
most of the modes contributing to the optimal growth belong to the vorticity
and entropy branches of the eigenspectrum. In fact, the phase speed
$c_{ph,opt}$, obtained as ratio between $\omega_{opt}$ and $\alpha_{opt}$, is
nearly 1. With respect to the value of $G_{opt}$ for both frameworks, a
difference of about $15\%$ is found. In conclusion, the same physical
mechanisms can be observed for both temporal and spatial analysis.

##  Appendix E The choice of energy norm

All non-modal optimizations presented in Sec. IV are performed with the new
energy norm introduced in Sec. II.3.2. This norm considers both the kinetic
and internal energy of a perturbation in a non-ideal gas flow. One can now
question whether the energy amplification is dependent on the choice of the
norm. To answer this question and analyze the sensitivity of the results in
Sec. IV, the energy norm in Eq. (21) is modified. When choosing only the
kinetic energy, i.e.,
$\mathbf{M}=\operatorname{diag}\left(0,\bar{\rho},\bar{\rho},\bar{\rho},0\right)$,
an infinite, and thus unphysical, energy amplification is achieved. Instead,
the internal energy content is intentionally altered by considering: (1)
frozen internal energy, i.e.,
$\mathbf{M}=\operatorname{diag}\left(1,\bar{\rho},\bar{\rho},\bar{\rho},1\right)$,
(2) internal energy according to ideal gas, e.g.,
$\mathbf{M}=\operatorname{diag}\left(R_{g}\bar{T}/\bar{\rho},\bar{\rho},\bar{\rho},\bar{\rho},\bar{\rho}\,c_{v}/(Ec_{\infty}\bar{T})\right)$.
Another alternative, hereafter (3), is to take the non-ideal optimal
perturbations obtained in Sec. II.3.3 and re-scale them by the kinetic energy
norm. Thus, following Eq. (27), the new optimal amplification for (3) is

$G_{(3)}=\max\dfrac{\|\mathbf{F}_{(3)}\bm{\Lambda}\bm{\kappa}(0)\|^{2}_{2}}{\|\mathbf{F}\bm{\kappa}(0)\|^{2}_{2}},\quad\bm{\kappa}(0)=\mathbf{F}^{-1}\mathbf{r},$
(39)

where $\mathbf{F}$, $\bm{\Lambda}$, $\mathbf{F}^{-1}$, $\mathbf{r}$ (right
singular eigenvector) are obtained with $\mathbf{M}$ as in Eq. (24), whereas
$\mathbf{F}_{(3)}$ is the Cholesky decomposition of
$\mathbf{A}_{(3)}=\mathbf{F}_{(3)}^{H}\mathbf{F}_{(3)}$ (with
$\mathbf{A}_{(3),kl}=\int_{0}^{\infty}\hat{\mathbf{q}}^{H}_{k}\mathbf{M}_{(3)}\hat{\mathbf{q}}_{l}\
dy$) as a function of the kinetic energy norm matrix
$\mathbf{M}_{(3)}=\operatorname{diag}\left(0,\bar{\rho},\bar{\rho},\bar{\rho},0\right)$.
Energy norms based on (1) and (2) were previously used in Ref. [7].
In order to account for the largest non-ideal effects on the norm, the two
transcritical cases of Tab. 1, with their optimal growth reported in Tab. 2,
are examined.
In Fig. 32(a), the maximum energy amplification $G_{max}(\omega)$ at constant
$\beta=\beta_{opt}$, i.e., $\beta_{opt}=0.45$, is shown for case T09w105. The
optimal energy amplification is found, independently of the energy norm
adopted, around $\omega=0.013$ (marked by a colored star ($\filledstar$)
symbol). The non-ideal energy norm yields the largest transient growth (blue
dash-dotted line). The shape of optimal perturbations displayed in Fig. 6 is
not affected by the energy-norm modifications (not presented here for the sake
of brevity). In Fig. 32(b), the maximum energy amplification $G_{max}(\beta)$
at constant $\omega=\omega_{opt}$, i.e., $\omega_{opt}=0$, is presented for
case T11w095. Here, we notice a greater dependence on the energy norm,
especially when the internal energy is kept constant (case (1), yellow
continuous line). In this case, there is an actual increase in the internal
energy compared to the other norms. Moreover, a shift of the maximum energy
amplification is notable. Similar to the transcritical T09w105 case, the shape
of the optimal perturbations remains intact except for the density
perturbations, where a smaller $\hat{\rho}$-amplitude in case (1) is
witnessed. Overall, energy amplifications with the re-scaled energy norm (case
(3), purple dotted line) show an analogous transient growth as the ideal-gas
norm cases (case (2), red dashed line).

Figure 32: Effect of the energy norm on the energy amplification: non-ideal
energy norm, ideal-gas energy norm (2), frozen internal energy norm (1), re-
scaled non-ideal energy norm by an kinetic energy norm (3). (a) T09w105, as a
function of $\omega$ for $\beta=\beta_{opt}$; (b) T11w095, as a function of
$\beta$ for $\omega=0$. The maximum energy amplification is indicated with a
colored star ($\filledstar$) symbol.

##  Appendix F Inviscid vorticity perturbation equation

The equation for the wall-normal vorticity perturbation $\omega^{\prime}_{y}$
is obtained hereafter. The wall-normal vorticity is split as
$\omega_{y}=\bar{\omega}_{y}+\omega^{\prime}_{y}$ and inserted into Eq. (29).
Subsequently, the following assumptions, in agreement with Sec. II.3, are
drawn: (1) 2-D base flow and locally parallel flow, i.e.,
$\bar{\omega}_{y}=\bar{v}=\bar{w}=0,\partial(\bar{\cdot})/\partial
x=0,\partial(\bar{\cdot})/\partial z=0$, (2) base-flow wall-normal pressure
gradient $\partial\bar{p}/\partial y$ is 0, (3) linearization, i.e.,
$q^{2}_{i}=q_{i}q_{j}=0$. After subtraction of the base flow, the baroclinic
term in Eq. (29) becomes:

$\dfrac{1}{\rho^{2}}\left(\dfrac{\partial\rho}{\partial z}\dfrac{\partial
p}{\partial x}-\dfrac{\partial\rho}{\partial x}\dfrac{\partial p}{\partial
z}\right)\quad\Longrightarrow\quad\dfrac{1}{\bar{\rho}^{2}}\left(\dfrac{\partial\rho^{\prime}}{\partial
z}\dfrac{\partial\bar{p}}{\partial x}+\dfrac{\partial\bar{\rho}}{\partial
z}\dfrac{\partial p^{\prime}}{\partial x}-\dfrac{\partial\bar{\rho}}{\partial
x}\dfrac{\partial p^{\prime}}{\partial
z}-\dfrac{\partial\rho^{\prime}}{\partial x}\dfrac{\partial\bar{p}}{\partial
z}\right)=0.$ (40)

Thus, the final expression of the wall-normal vorticity perturbation equation
is

$\dfrac{\partial\omega^{\prime}_{y}}{\partial
t}+\bar{u}\dfrac{\partial\omega^{\prime}_{y}}{\partial x}=-\dfrac{\partial
v^{\prime}}{\partial z}\dfrac{\partial\bar{u}}{\partial y},$ (41)

with absent baroclinic influence.

##  Appendix G Modal stability analysis: oblique perturbations

The influence of 3-D perturbations on modal growth needs to be investigated
for the $N$-factor comparison in Sec. IV.4. In Ref. [8], the maximum local
growth ratio, defined as
$\max_{\forall\omega}\\{\alpha_{i,\beta=\beta_{max}}/\alpha_{i,\beta=0}\\}$,
was calculated for an adiabatic wall at finite Eckert numbers. 2-D
perturbations were found to be most unstable in the subcritical regime,
whereas a maximum local-growth ratio larger than one was obtained for oblique
disturbances in the supercritical regime. In the transcritical regime, Mode I
and Mode II were most locally amplified in 3-D and 2-D, respectively.
In this study, a similar investigation is undertaken with a special focus on
the integral amplification rather than the local one. Stability diagrams are
first obtained in the $Re_{\delta}$–$\omega$–$\beta$ space up to
$Re_{\delta}=2000$. Cuts of $Re_{\delta}$–$\omega$ are then showcased at
different $\beta$ with iso-contours of the largest $N$-factor over all
spanwise wavenumbers. Figs. 33(a,c) and Figs. 33(b,d) exemplify non-
transcritical cases T09w085 and T11w105, respectively. Their largest growth
rate in the stability diagram and their largest integral amplification are
detected at $\beta=\beta_{max}=0$, with larger values for the subcritical
wall-cooling regime. For instance, the largest eigenvalues correspond to
$(c_{r},\alpha_{i},Re_{\delta},\beta)=(0.301,-0.010,1400,0)$ and
$(c_{r},\alpha_{i},Re_{\delta},\beta)=(0.252,-0.0053,2200,0)$ for T09w085 and
T110w105, respectively. The real part of the phase velocity $c_{r}$ is
computed as $\omega/\alpha_{r}$. Note that, as illustrated in Figs. 33(c,d),
there are spanwise wavenumbers at which 3-D modes are locally more amplified
than 2-D modes. Differently than at finite Eckert numbers [8], all non-
transcritical cases at $M_{\infty}=10^{-3}$ have their maximum local-growth
ratio and largest integral amplification for 2-D perturbations.

Figure 33: Slices of $-\alpha_{i}$-contours in the
$Re_{\delta}$–$\omega$–$\beta$ (a,b) and neutral stability curves in the
$Re_{\delta}$–$\omega$ (c,d) space: (a,c) T09w085, (b,d) T11w105. Arrows in
(c,d) indicate a $\beta$-increase ($\Delta\beta=0.05$) of branch I and II
(from solid blue to solid black line). Maximum $N$-factors are indicated in
blue and are located always at $\beta=\beta_{max}=0$.

The effect of an increasing spanwise wavenumber is displayed in Fig. 34 for
the transcritical case T09w105. Mode II exhibits its largest growth rate over
all frequencies at $\beta=\beta_{m}=0$, while Mode I can also be most unstable
for 3-D disturbances. The largest eigenvalues are found at
$(c_{r},\alpha_{i},Re_{\delta},\beta)=(0.457,-0.011,2000,0)$ and
$(c_{r},\alpha_{i},Re_{\delta},\beta)=(0.312,-0.072,2060,0)$ for Mode I and
Mode II, respectively. Nevertheless, Mode II is almost unaffected by a
$\beta$-increase, being extremely unstable in the considered parameter space.
Hence, $N$-factors up to $90$ are achieved at $Re_{\delta}=2000$. For Mode I,
an increase in $\beta$ leads to a stabilization analogous to the non-
transcritical cases.

Figure 34: Case T09w105. Slices of $-\alpha_{i}$-contours in the
$Re_{\delta}$–$\omega$–$\beta$ (a,b) and neutral stability curves in the
$Re_{\delta}$–$\omega$ (c,d) space: (a,c) Mode I, (b,d) Mode II. Arrows in
(c,d) indicate a $\beta$-increase ($\Delta\beta=0.05$) of branch I and II
(from solid blue to dashed line in (c) and to dash-dotted line in (d)).
Maximum $N$-factors are indicated in blue and are located always at
$\beta=\beta_{max}=0$.

For case T11w095, in Fig. 35, only one highly unstable mode is detected in
agreement with Ref. [12]. It is inviscid, given the GIP in Fig. 2(c). Similar
to Mode II in the transcritical wall-heating case, this mode is barely
unaffected by a $\beta$-increase. Yet, it extends up to the flat-plate leading
edge and to relatively high frequencies. Interestingly, the largest growth
rate, i.e., $\alpha_{i,max}=-0.0946$, is found at a very low local Reynolds
number of $Re_{\delta}=220$ with a phase speed of $c_{r}=0.217$. With respect
to the integral amplification, we observe a $N$-factor of 15 already at
$Re_{\delta}=200$.
We can conclude that, for all cases in Tab. 1, the current modal analysis has
revealed the ”dominance” of 2-D disturbances at $M=10^{-3}$. Both the local
amplification (over most of the frequency spectrum) as well as the maximum
local-growth ratio, and the $N$-factor are found largest for 2-D modes. Thus,
these results clearly indicate that reducing the Mach number has the effect of
shifting the maximum (local and integral) amplification from a 3-D to a 2-D
mode independently of the considered thermodynamic regime at supercritical
pressure.

Figure 35: Case T11w095. Slices of $-\alpha_{i}$-contours in the
$Re_{\delta}$–$\omega$–$\beta$ (a) and neutral stability curves in the
$Re_{\delta}$–$\omega$ (b) space. Arrows in (b) indicate a $\beta$-increase
($\Delta\beta=0.05$) of branch I (from solid blue to dash-dotted line).
Maximum $N$-factors are indicated in blue and are located always at
$\beta=\beta_{max}=0$.

## References

  * Liao  _et al._ [2019] G. Liao, L. Liu, E. Jiaqiang, F. Zhang, J. Chen, Y. Deng, and H. Zhu, Effects of technical progress on performance and application of supercritical carbon dioxide power cycle: a review, Energy Convers. Manag. 199, 111986 (2019).
  * Kawai [2019] S. Kawai, Heated transcritical and unheated non-transcritical turbulent boundary layers at supercritical pressures, J. Fluid Mech. 865, 563–601 (2019).
  * Peeters  _et al._ [2016] J. W. R. Peeters, R. Pecnik, M. Rohde, T. H. J. J. van der Hagen, and B. J. Boersma, Turbulence attenuation in simultaneously heated and cooled annular flows at supercritical pressure, J. Fluid Mech. 799, 505–540 (2016).
  * Yoo [2013] J. Y. Yoo, The turbulent flows of supercritical fluids with heat transfer, Annu. Rev. Fluid Mech. 45, 495–525 (2013).
  * Nemati  _et al._ [2016] H. Nemati, A. Patel, B. J. Boersma, and R. Pecnik, The effect of thermal boundary conditions on forced convection heat transfer to fluids at supercritical pressure, J. Fluid Mech. 800, 531–556 (2016).
  * Gloerfelt  _et al._ [2020] X. Gloerfelt, J.-C. Robinet, L. Sciacovelli, P. Cinnella, and F. Grasso, Dense-gas effects on compressible boundary-layer stability, J. Fluid Mech. 893, A19 (2020).
  * Ren  _et al._ [2019a] J. Ren, S. Fu, and R. Pecnik, Linear instability of Poiseuille flows with highly non-ideal fluid, J. Fluid Mech. 859, 89–125 (2019a).
  * Ren  _et al._ [2019b] J. Ren, O. Marxen, and R. Pecnik, Boundary-layer stability of supercritical fluids in the vicinity of the Widom line, J. Fluid Mech. 871, 831–864 (2019b).
  * Ly and Ihme [2022] N. Ly and M. Ihme, Destabilization of binary mixing layer in supercritical conditions, J. Fluid Mech. 945, R2 (2022).
  * Bugeat  _et al._ [2022a] B. Bugeat, P. C. Boldini, and R. Pecnik, On the new unstable mode in the boundary layer flow of supercritical fluids, in  _Proceedings of the 12th International Symposium on Turbulence and Shear Flow Phenomena (TSFP-12)_ (2022).
  * Bugeat  _et al._ [2022b] B. Bugeat, P. C. Boldini, A. M. Hasan, and R. Pecnik, Instability in strongly stratified plane Couette flow, with application to supercritical fluids, arXiv preprint arXiv:2211.04935 (2022b).
  * Ren and Kloker [2022] J. Ren and M. Kloker, Instabilities in three-dimensional boundary-layer flows with a highly non-ideal fluid, J. Fluid Mech. 951, A9 (2022).
  * Banuti [2015] D. T. Banuti, Crossing the Widom-line – supercritical pseudo-boiling, J. Supercrit. Fluids 98, 12–16 (2015).
  * Schmid [2007] P. J. Schmid, Nonmodal stability theory, Annu. Rev. Fluid Mech. 39, 129 (2007).
  * Reshotko [2001] E. Reshotko, Transient growth: A factor in bypass transition, Phys. Fluids 5, 1067–1075 (2001).
  * Reshotko and Tumin [2000] E. Reshotko and A. Tumin, The blunt body paradox – a case for transient growth, in  _Laminar-Turbulent Transition_ (Springer, Berlin, 2000) p. 403–408.
  * Hanifi  _et al._ [1996] A. Hanifi, P. J. Schmid, and D. S. Henningson, Transient growth in compressible boundary layer flow, Phys. Fluids 8, 826–837 (1996).
  * Andersson  _et al._ [1999] P. Andersson, M. Berggren, and D. S. Henningson, Optimal disturbances and bypass transition in boundary layers, J. Fluid Mech. 11, 134–150 (1999).
  * Landahl [1980] M. T. Landahl, A note on an algebraic instability of inviscid parallel shear flows, J. Fluid Mech. 98, 243–251 (1980).
  * Luchini [2000] P. Luchini, Reynolds number independent instability of the boundary layer over a flat surface: optimal perturbations, J. Fluid Mech. 404, 289–309 (2000).
  * Hanifi and Henningson [1998] A. Hanifi and D. S. Henningson, The compressible inviscid algebraic instability for streamwise independent disturbances, Phys. Fluids 10 (1998).
  * Tumin and Reshotko [2001] A. Tumin and E. Reshotko, Spatial theory of optimal disturbances in boundary layers, Phys. Fluids 13, 2097–2104 (2001).
  * Bitter [2015] N. P. Bitter,  _Stability of hypervelocity boundary layers_ , Ph.D. thesis, California Institute of Technology (2015).
  * Paredes  _et al._ [2016] P. Paredes, M. M. Choudhari, F. Li, and C. L. Chang, Optimal growth in hypersonic boundary layers, AIAA J. 54, 3050–3061 (2016).
  * Govindarajan and Sahu [2014] R. Govindarajan and K. C. Sahu, Instabilities in viscosity-stratified flow, Ann. Rev. Fluid Mech. 46, 331–353 (2014).
  * Malik  _et al._ [2008] M. Malik, J. Dey, and M. Alam, Linear stability, transient energy growth, and the role of viscosity stratification in compressible plane Couette flow, Phys. Rev. E 77, 036322 (2008).
  * Saikia  _et al._ [2017] B. Saikia, A. Ramachandran, K. Sinha, and R. Govindarajan, Effects of viscosity and conductivity stratification on the linear stability and transient growth within compressible Couette flow, Phys. Fluids 29, 1–20 (2017).
  * Chikkadi  _et al._ [2005] V. Chikkadi, A. Sameen, and R. Govindarajan, Preventing transition to turbulence: a viscosity stratification does not always help, Phys. Rev. Lett. 95, 264504 (2005).
  * Sameen and Govindarajan [2007] A. Sameen and R. Govindarajan, The effect of wall heating on instability of channel flow, J. Fluid Mech. 577, 417–442 (2007).
  * Sameen  _et al._ [2011] A. Sameen, R. Bale, and R. Govindarajan, The effect of wall heating on instability of channel flow – CORRIGENDUM, J. Fluid Mech. 673, 603–605 (2011).
  * Jose  _et al._ [2020] S. Jose, L. Brandt, and R. Govindarajan, Localisation of optimal perturbations in variable viscosity channel flow, Int. J. Heat Fluid Flow 85, 108588 (2020).
  * Parente  _et al._ [2020] E. Parente, J. C. Robinet, P. D. Palma, and S. Cherubini, Modal and nonmodal stability of a stably stratified boundary layer flow, Phys. Rev. Fluids 5, 113901 (2020).
  * Robinet and Gloerfelt [2019] J.-C. Robinet and X. Gloerfelt, Instabilities in non-ideal fluids, J. Fluid Mech. 880, 1–4 (2019).
  * Lemmon  _et al._ [2013] E. W. Lemmon, M. L. Huber, and M. O. Mclinden, NIST Standard Reference Database 23: Reference Fluid Thermodynamic and Transport Properties - REFPROP, Version 9.1 (2013), Available at: http://www.nist.gov/srd/nist23.cfm.
  * White [2006] F. M. White,  _Viscous Fluid Flow, 3rd Edition_ (McGraw-Hill, Boston, 2006).
  * Chu [1965] T. B. Chu, On the energy transfer to small disturbances in fluid flow (part I), Acta Mech. 1, 215–234 (1965).
  * Mack [1969] L. M. Mack, Boundary layer stability theory, JPL Rep. 900-227. Jet Propulsion Laboratory. (1969).
  * Chen  _et al._ [2022] X. Chen, L. Wang, and S. Fu, Energy transfer of hypersonic and high-enthalpy boundary layer instabilities and transition, Phys. Rev. Fluids 7, 033901 (2022).
  * Schmid and Henningson [1994] P. J. Schmid and D. S. Henningson, Optimal energy density growth in Hagen–Poiseuille flow, J. Fluid Mech. 277, 197–225 (1994).
  * Corbett and Bottaro [2000] P. Corbett and A. Bottaro, Optimal perturbations for boundary layers subject to stream-wise pressure gradient, Phys. Fluids 12, 120–130 (2000).
  * Butler and Farrell [1992] K. M. Butler and B. F. Farrell, Three-dimensional optimal perturbations in viscous shear flow, Phys. Fluids A 4, 1637–1650 (1992).
  * Brandt [2014] L. Brandt, The lift-up effect: the linear mechanism behind transition and turbulence in shear flows, Eur. J. Mech. (B/Fluids) 47, 80–96 (2014).
  * Farrell and Ioannou [1993] B. F. Farrell and P. J. Ioannou, Optimal excitation of three-dimensional perturbations in viscous constant shear flow, Phys. Fluids A 5, 1390–1400 (1993).
  * Hack and Moin [2017] M. J. P. Hack and P. Moin, Algebraic disturbance growth by interaction of Orr and lift-up mechanisms, J. Fluid Mech. 829, 112–126 (2017).
  * Orr [1907] W. M. F. Orr, The stability or instability of the steady motions of a perfect liquid and of a viscous liquid, Proc. R. Irish Acad. A 27, 9–138 (1907).
  * Berlin  _et al._ [1994] S. Berlin, A. Lundbladh, and D. Henningson, Spatial simulations of oblique transition in a boundary layer, Phys. Fluids 6, 1949–1951 (1994).
  * Tempelmann  _et al._ [2012] D. Tempelmann, A. Hanifi, and D. S. Henningson, Spatial optimal growth in three-dimensional compressible boundary layers, J. Fluid Mech. 704, 251–279 (2012).
  * Gustavsson [1991] L. H. Gustavsson, Energy growth of three-dimensional disturbances in plane Poiseuille flow, J. Fluid Mech. 224, 241–260 (1991).
  * Jiao  _et al._ [2021] Y. Jiao, Y. Hwang, and S. I. Chernyshenko, Orr mechanism in transition of parallel shear flow, Phys. Rev. Fluids 6, 023902 (2021).
  * Tumin [2007] A. Tumin, Three-dimensional spatial normal modes in compressible boundary layers, J. Fluid Mech. 586, 295–322 (2007).
  * Kaminski  _et al._ [2014] A. K. Kaminski, C. P. Caulfield, and J. R. Taylor, Transient growth in strongly stratified shear layers, J. Fluid Mech. 758, R4 (2014).
  * Tumin and Reshotko [2003] A. Tumin and E. Reshotko, Optimal disturbances in compressible boundary layers, AIAA J. 41, 2357–2363 (2003).
  * Levin and Henningson [2003] O. Levin and D. S. Henningson, Exponential vs algebraic growth and transition prediction in boundary layer flow, Flow Turbul. Combust. 70, 183–210 (2003).
  * Reshotko and Tumin [2004] E. Reshotko and A. Tumin, Role of transient growth in roughness-induced transition, AIAA J. 42 (4), 766–770 (2004).


#  Robustness Analysis of Classification Using Recurrent Neural Networks with
Perturbed Sequential Input

Guangyi Liu, Arash Amini, Martin Takáč, and Nader Motee G.L., A.A., and N.M.
are with the Department of Mechanical Engineering and Mechanics, Lehigh
University, Bethlehem, PA 18015, USA<EMAIL_ADDRESS>is
with the Machine Learning Department at Mohamed bin Zayed University of
Artificial Intelligence (MBZUAI), Masdar City, Abu Dhabi, United Arab Emirates
<EMAIL_ADDRESS>
###### Abstract

For a given stable recurrent neural network (RNN) that is trained to perform a
classification task using sequential inputs, we quantify explicit robustness
bounds as a function of trainable weight matrices. The sequential inputs can
be perturbed in various ways, e.g., streaming images can be deformed due to
robot motion or imperfect camera lens. Using the notion of the Voronoi diagram
and Lipschitz properties of stable RNNs, we provide a thorough analysis and
characterize the maximum allowable perturbations while guaranteeing the full
accuracy of the classification task. We illustrate and validate our
theoretical results using a map dataset with clouds as well as the MNIST
dataset.

##  I Introduction

Real-time perception and classification have been among the most exciting
topics in computer vision-based and machine learning-based robotic
applications. However, in the most perception-based applications, the
uncertainties and perturbations prevail in every section of a learned model
from the input perturbation [6, 24] to the numerical error [28]. In order to
ensure reliable performance for these applications, some inevitable questions
need to be answered: (i) If input perturbation exists, does the learned
framework still exhibit robust performance? (ii) What conditions does the
learned model needs to satisfy to achieve robustness? (iii) Does there exist
characteristics that can measure the robustness of the learned model? Finding
the answer to these questions will significantly facilitate tackling
perception-based problems since most learned models suffer from the fragility
of the input perturbations [30].

In many applications involving area coverage, consensus, map classification,
rendezvous [2], the robot can only sample the localized information of the
environment [29], i.e., partially observation, at each time step. In order to
acquire adequate observations, the robot could consider traversing the
environment while collecting local information as sequential data and
attempting to learn their inter-correlation with recurrent neural networks.
The same type of framework that uses localized observations and recurrent
neural networks to learn the image and map classification is illustrated in
our previous works [22, 21, 16, 17], and shows promising performance in
various scenarios.

In this paper, we consider the map (image) classification problem with the
sampled sequential images as a motivational example, which is illustrated in
Fig. 1. Instead of treating the learned classification model as a black box
and solely focusing on its performance, we consider how it will perform with
perturbed input and analyze its robustness in terms of its learned weights
from the neural networks. The robustness analysis is achieved by considering
the stability properties of the recurrent models that learn the
interconnection of past observations and quantifying the classification
criterion via the Voronoi partitioning to obtain the robustness conditions.

Figure 1: The above figure depicts the motivational application of map
classification with the aerial robot. The robot traverses the environment and
collects the visual inputs as a sequence (order denoted by arrows).

The robustness of a classification model can be considered as not miss-
classifying under the adversary attacks or perturbations [18], and recent
research has made progress on investigating and ensuring the robustness of
neural network models under the adversarial attack [31, 5]. Our work exhibits
the novelty and differences to these works as: Instead of considering the
Boolean classifier, we propose a robustness analysis of multi-class classifier
[19], in which a novel representation with Voronoi diagram [4] is used to
construct a quantifiable classification criterion. To ensure the robustness of
a learned model, we implement the similar ideas of using expert demonstrations
[25], for which we require the perturbed result to stay close to the nominal
results, see §VI. To evaluate the deviation caused by the perturbation at the
output stage, we seek the boundedness for every section of the classification
model. For instance, inspired by the recent work on constructing a bound for
neural network models, to name a few: [20, 12, 31], we seek similar
boundedness for RNNs, which will provide an error estimate when the statistics
of the input perturbation is available [1].

Our Contributions:  In this paper, we construct a formal approach to analyze
the robustness of the recurrent multi-label classification framework with
localized sequential inputs when there exist input perturbations. Furthermore,
we also quantify the multi-label classification criterion that uses the
$\arg\max$ function. Our analysis shows that a robust classification with full
accuracy is guaranteed when the RNNs are stable and the input perturbation is
below the maximum allowable deviation. These results motivate us to turn our
research efforts to explore further how the classification robustness can show
its effect in a closed-loop model where the robot can select its sampling
routine based on its observations.

The rest of the paper is organized as follows. In §III, we introduce the
problem setting and the classification model. The possible origins of
perturbations are illustrated in §IV. The error estimates and the boundedness
of stable RNNs are presented in §V. Our main result is presented in §VI, where
the quantifiable classification criterion and robust classification conditions
are presented. The theoretical findings are validated in §VII by simulations
in both MNIST dataset [14] and Campus Map dataset [17].

##  II Mathematical Notations

The $n-$dimensional Euclidean space with elements
$\bm{z}=[z_{1},\dots,z_{n}]^{T}$ is denoted by $\mathbb{R}^{n}$, where
$\mathbb{R}_{+}$ will denote the positive orthant of $\mathbb{R}^{n}$. The set
of standard Euclidean basis for $\mathbb{R}^{n}$ is represented by
$\\{\bm{e}_{1},\dots,\bm{e}_{n}\\}$. We denote the $n\times n$ identity matrix
as $I$ and the vector of all ones as $\bm{1}_{n}$, respectively. The $i$’th
element of a vector $x$ is shown by $x_{i}$ and the $i$’th row of a matrix $A$
is represented by $(A)_{i}$. The induced matrix norm by vector norm $\|\hskip
1.42271pt.\hskip 1.42271pt\|$ is also shown by $\|A\|$ [26]. Let us define the
collection of all feasible probability vectors $p$ [15] as
$\mathcal{P}_{m}=\\{p\in\mathbb{R}_{+}^{m}~{}|~{}p^{T}\bm{1}_{m}=1\\}$. For a
sequence of vectors $x=(x(1),x(2),\ldots,x(T))$, the $\ell_{\infty}$-norm of
$x$ is defined by

$\|x\|_{\ell_{\infty}}=\max_{t\in\\{1,\dots,T\\}}\|x(t)\|.$ (1)

##  III Problem Statement

Figure 2: This figure shows the image (map) classification model at time $t$.

Suppose there exist $m$ unique pre-labeled environments for classification
purposes. A robot is deployed into the environment to collect samples for
classification. The sample, e.g., a vector that contains multiple states of
the environment, from data sequence $x$. The data sequence $x$ with length $T$
can be represented in terms of its components as $x=(x(1),x(2),\dots,x(T))$.
Let us denote by $\mathcal{X}_{k}^{\prime}$ the set which contains adequate
sampled data sequences $x$ from the $k$’th environment, i.e., the training
set. Our proposed classifier, depicted in Fig. 2, is trained with
$\mathcal{X}_{k}^{\prime}$ for all $k=1,...,m$ to classify the label of $x$,
which is sampled from an environment with unknown labels. In the classifier,
each individual sample of $x$ are fed recurrently to a RNN model [1], whose
dynamics can be represented in a compact form by

$\displaystyle h(t;x)=F\left(h(t-1;x),x(t)\right),$ (2)

where $h(t;x)\in\mathbb{R}^{b}$ is the state and $x(t)\in\mathbb{R}^{a}$ is
the input. The (history) state $h(t;x)$ memorizes all the past information
about inputs $x$ up to time $t$. The terminal state $h(T;x)$ is used for
classification by passing it through a fully-connected layer

$\displaystyle q(T;x)=W_{c}\>h(T;x)+b_{c}$ (3)

and a Softmax function

$\displaystyle p(T;x)=\texttt{Softmax}\left(q(T;x)\right),$ (4)

where the weight matrix $W_{c}\in\mathbb{R}^{m\times b}$ and the bias vector
$b_{c}\in\mathbb{R}^{m}$ are trainable. The belief vector $p(T;x)$ is utilized
to represent the classification result from the sequence $x$.

Our objective is to provide a thorough robustness analysis of the
classification model using stable RNNs and quantify their robustness bounds in
terms of their trainable weight matrices.

##  IV Origins of Perturbations

Figure 3: An input can be perturbed in various ways. This figures illustrate
some possible cases for the map classification task.

There exist several ways by which data can be collected as a sequence from the
environment; for example, an agent can navigate in an environment and take
localized observations [16] or a camera with a fixed location can capture
images from a time-varying scene [9]. In most real-world applications, such
raw observations are pre-processed, e.g., by neural networks, and the data
sequence already carries relevant features of the observed raw data. Let us
assume that this pre-processing can be modeled by a nonlinear map

$x(t)=\Phi\left(O(t)\right),$

where $O(t)\in\mathcal{O}$ denotes the raw observation sampled at time $t$,
and $\mathcal{O}$ is the space of all observables.

Perturbations can affect the data quality through several possible sources
during the sampling process. The observer may deviate from its initially
planned sampling routine due to dynamic noise in its motion planning [7] and
take (slightly) deviated samples from nearby scenes in the environment. The
raw observations may lose quality due to the environmental noise [3], e.g.,
change in light intensity, cloudiness, and blurriness. The raw observations
may also experience various types of deformation, e.g., camera rotation or
distortion [30]. Fig. 3 depicts some of these perturbations. The effect of
uncertainty in all these cases can be modeled by

$\tilde{x}(t)=\Phi\left(\tau\left(O(t)\right)+\xi(t)\right),$ (5)

where $\tilde{x}(t)$ represents the perturbed data,
$\xi:\mathbb{R}\rightarrow\mathcal{O}$ is an additive bounded stochastic noise
or a bounded deterministic disturbance, and
$\tau:\mathcal{O}\rightarrow\mathcal{O}$ is a deformation map that models
sample deformation due to sensor movements (e.g., translation and rotation,
and scaling).

##  V Stable Recurrent Neural Networks and Their Error Estimates

In order to analyze the robustness of the classification model, we will
evaluate how the input perturbation is affecting the classification result.
Let us first identify what input sequences can generate correct
classifications. Recall that the model is trained with sequences from
$\mathcal{X}_{k}^{\prime}$ for all $k=1,...,m$, and not every
$x\in\mathcal{X}_{k}^{\prime}$ will carry enough information to reveal the
environment (e.g., a sequence consists of repeating scenes). Hence, only a
subset of the sequences from the training set may generate the correct
classification. We represent those data sequences that generate correct
classifications by $\mathcal{X}_{k}^{*}$, where
$\mathcal{X}_{k}^{*}\subseteq\mathcal{X}_{k}^{\prime}$. This relation is
depicted in Fig. 4.

To measure how the input perturbation will affect the classification, let us
consider a nominal sequence $x\in\mathcal{X}_{k}^{*}$ and its corresponding
perturbed sequence $\tilde{x}\in\mathcal{X}_{k}$ obtained from (5). The set
$\mathcal{X}_{k}$ is the space of all possible sequences for the $k$’th class.
In this section, we aim to evaluate the deviation generated by $x$ and
$\tilde{x}$ in terms of belief vectors, i.e., $p(T;x)$ and $p(T;\tilde{x})$.
The first step is to consider the deviation generated at the output of the
RNN, i.e., $h(t;x)$ and $h(t;\tilde{x})$.

In order to demonstrate our next result, let us introduce the concept of
stable RNNs and its Lipschitz property. A stable RNN [20] provides the
boundedness to its output when the input vectors are identical, i.e.,
$x(t)=\tilde{x}(t)$.

Figure 4: This figure depicts the relation of the space of all possible
sequences $\mathcal{X}_{k}$, the training set $\mathcal{X}_{k}^{\prime}$, and
the nominal set $\mathcal{X}_{k}^{*}$, i.e.,
$\mathcal{X}_{k}^{*}\subseteq\mathcal{X}_{k}^{\prime}\subseteq\mathcal{X}_{k}$.

######  Definition 1.

A recurrent neural network model (2) is stable (contractive) if there exists a
constant $\lambda\in(0,1)$ such that for any
$h(t-1;x),h(t-1;\tilde{x})\in\mathbb{R}^{b}$,

$\|h(t;x)-h(t;\tilde{x})\|\leq\lambda\,\|h(t-1;x)-h(t-1;\tilde{x})\|,$ (6)

where $h(t;x)=F\big{(}h(t-1;x),x(t)\big{)}$ and
$h(t;\tilde{x})=F\big{(}h(t-1;\tilde{x}),x(t)\big{)}$.

It is also known that the RNN model is Lipschitz continuous with respect to
its input [20, 12]. Then, for every two sequences $x$ and $\tilde{x}$ with
identical history state at $t-1$, i.e., $h(t-1;x)=h(t-1;\tilde{x})$, one has

$\displaystyle\|h(t;x)-h(t;\tilde{x})\|\leq\kappa~{}\|x(t)-\tilde{x}(t)\|,$
(7)

where $\kappa\in\mathbb{R}_{+}$ is the real Lipschitz constant,
$h(t;x)=F\big{(}h(t-1;x),x(t)\big{)}$, and
$h(t;\tilde{x})=F\big{(}h(t-1;x),\tilde{x}(t)\big{)}$. Then, the following
result shows the upper bound for the deviation between the nominal belief
$p(T;x)$ to the perturbed belief $p(T;\tilde{x})$ given a stable RNN model.

######  Theorem 1.

For the classification model with a stable RNN model, one has

$\|p(T;x)-p(T;\tilde{x})\|\leq\eta\,\|x-\tilde{x}\|_{\ell_{\infty}}$

for every two sequences $x\in\mathcal{X}_{k}^{*}$ and
$\tilde{x}\in\mathcal{X}_{k}$, where

$\eta=\frac{\kappa~{}\|W_{c}\|}{(1-\lambda)\sqrt{m}}.$ (8)

###### Proof.

The initial hidden states are set to $h(0;x)=0$ and $h(0;\tilde{x})=0$ in both
training and testing stage. Considering the Lipschitz continuity (7), the
history states generated by sequences $x$ and $\tilde{x}$ at $t=1$ follows

$\displaystyle\|h(1;x)-h(1;\tilde{x})\|\leq\kappa\|x(1)-\tilde{x}(1)\|.$ (9)

For the next time step $t=2$, we can obtain (10) by using the triangle
inequality

$\|h(2;x)-h(2;\tilde{x})\|=\|h(2;x)-F\big{(}h(1;x),\tilde{x}(2)\big{)}\\\
+F\big{(}h(1;x),\tilde{x}(2)\big{)}-h(2;\tilde{x})\|\\\
\leq\|F\big{(}h(1;x),x(2)\big{)}-F\big{(}h(1;x),\tilde{x}(2)\big{)}\|\\\
+\|F\big{(}h(1;x),\tilde{x}(2)\big{)}-F\big{(}h(1;\tilde{x}),\tilde{x}(2)\big{)}\|.$
(10)

The first half of (10) can be bounded using (7),

$\displaystyle\|F\big{(}h(1;x),x(2)\big{)}-F\big{(}h(1;x),\tilde{x}(2)\big{)}\|$
$\displaystyle\leq\kappa~{}\|x(2)-\tilde{x}(2)\|$
$\displaystyle\leq\kappa~{}\|x-\tilde{x}\|_{\ell_{\infty}}.$

The second half of (10) can be bounded using (6) and then (9),

$\displaystyle\|F(h(1;x),\tilde{x}(2))-F(h(1;\tilde{x}),\tilde{x}(2))\|$
$\displaystyle\hskip 56.9055pt\leq\lambda~{}\|h(1;x)-h(1;\tilde{x})\|$
$\displaystyle\hskip
56.9055pt\leq\lambda\kappa~{}\|x(1)-\tilde{x}(1)\|\leq\lambda\kappa~{}\|x-\tilde{x}\|_{\ell_{\infty}}.$

Summarizing the above inequalities, the deviation of history states at $t=2$
obtains the following upper bound,

$\displaystyle\|h(2;x)-h(2;\tilde{x})\|\leq(1+\lambda)\,\kappa\,\|x-\tilde{x}\|_{\ell_{\infty}}.$

Repeating the above steps up to $t$, the deviation of
$\|h(t;x)-h(t;\tilde{x})\|$ is upper bounded by

$\displaystyle\|h(t;x)-h(t;\tilde{x})\|$
$\displaystyle\leq(1+\lambda+...+\lambda^{t-1})\kappa~{}\|x-\tilde{x}\|_{\ell_{\infty}}$
$\displaystyle\leq\frac{1-\lambda^{t}}{1-\lambda}\kappa~{}\|x-\tilde{x}\|_{\ell_{\infty}}.$

Given that $0<\lambda<1$, the terminal deviation at $T$ obtains the following
boundedness,

$\displaystyle\|h(T;x)-h(T;\tilde{x})\|\leq\frac{\kappa}{1-\lambda}~{}\|x-\tilde{x}\|_{\ell_{\infty}}.$
(11)

The similar bounds can be obtained for the fully connected layer (3),

$\displaystyle\|q(T;x)-q(T;\tilde{x})\|\leq\|W_{c}\|\frac{\kappa}{1-\lambda}~{}\|x-\tilde{x}\|_{\ell_{\infty}}.$

Then, the results follows immediately by considering the fact that the Softmax
function is $1/\sqrt{m}$ Lipshcitz continuous [13] with respect to its input
$q\in\mathbb{R}^{m}$. ∎

The above theorem asserts that, in the classifier, if the RNN model is stable,
the terminal belief difference is bounded by the maximum deviation along the
input sequences. This result provides the knowledge of the classifier’s
robustness by identifying an upper bound for the error generated by the
perturbed and the nominal data.

##  VI Classification Criterion and Convergence Conditions

Figure 5: This figure depicts the Voronoi partition of the classification sets
$\mathcal{P}_{m}$ for both $m=2$ (left) and $m=3$ (right).

In this part, we aim to represent the $\arg\max$ classification criterion in a
quantifiable manner and use the previously obtained deviation bounds in belief
vectors to perform the robustness analysis.

The classification process with $x$ is accomplished by identifying $\arg\max
p(T;x)$ as the class label, i.e., it will conclude as the $k$’th class if and
only if $p(T;x)_{k}>p(T;x)_{j}$ for all $j=1,\dots,m$ and $j\neq k$. This
criterion can be explicitly represented via the Voronoi partitioning [11]. Let
us denote by $V_{1},V_{2},\dots,V_{m}$ the Voronoi partition of the
probability vector space $\mathcal{P}_{m}$, i.e.,

$V_{k}=\left\\{p\in\mathcal{P}_{m}~{}|~{}p_{k}>p_{j},~{}\forall j\neq
k\right\\}.$ (12)

Some examples of the Voronoi partitioning is shown in Fig. 5. The above
partition is equivalent to the $\arg\max$ classification criterion, i.e.,
$\arg\max p(T;x)=p(T;x)_{k}$ if and only if $p(T;x)\in V_{k}$.

The next step is to convert the deviation $\|p(T;x)-p(T;\tilde{x})\|$ into a
Boolean classification result, i.e., “true” or “false”. To reveal our next
result, let us introduce the distance function between vectors and sets.

######  Definition 2.

The distance between two vectors $p,p^{\prime}\in\mathcal{P}_{m}$ is defined
as

$d(p,p^{\prime})=\|p-p^{\prime}\|$

and the distance between a vector and a set is defined as

$d(p,V_{j})=\inf_{p^{\prime}\in V_{j}}d(p,p^{\prime}).$

Let us express the $\arg\max$ classification criterion equivalently using the
distance function.

######  Lemma 1.

A belief vector $p\in\mathcal{P}_{m}$ will be classified as the $k$’th class
if and only if

$d(p,V_{k})=0,$ (13)

or if and only if

$d(p,V_{j})>0\text{ for all }j\neq k.$ (14)

###### Proof.

Suppose there exists a $p^{\prime}\in V_{k}$ such that $d(p,p^{\prime})=0$,
then one has $p=p^{\prime}$ and $p\in V_{k}$. On the other hand, if $p\in
V_{k}$, then we have $d(p,V_{k})\leq d(p,p)=0$ and $d(p,V_{k})=0$ by
definition.

For the equivalence, since $V_{k}\bigcap V_{j}=\emptyset$ for all $j\neq k$,
one has if $d(p,V_{k})=0$, then $d(p,V_{j})>0$ for any $V_{j}$. On the other
hand, consider the fact that

$V_{k}\subset\mathcal{P}_{m}\setminus\bigcup_{j\neq k}V_{j},$

one has if $d(p,V_{j})>0$, then $p\notin V_{j}$ and $p\in V_{k}$ 111The cases
of $p$ located on the boundary of the Voronoi partition is considered trivial
since its probability of happening is $0$.. ∎

The above result introduces a quantifiable classification criterion, which can
be used in all classification problems using the $\arg\max$ classifier. Given
a belief vector $p\in\mathcal{P}_{m}$, one can obtain the classification
result by checking if (13) or (14) are satisfied.

Using Lemma 1, we can establosh a connection between
$\|p(T;x)-p(T;\tilde{x})\|$ and the robustness of the classification model.
For a robust classification model, we expect $p(T;\tilde{x})$ to be located
within $V_{k}$. This relation can be validated by comparing $d(p(T;x),V_{j})$
with $\|p(T;x)-p(T;\tilde{x})\|$ for all $j=1,...,m$ and $j\neq k$. In order
to accomplish the analysis for an arbitrary perturbed sequence $\tilde{x}$,
let us introduce the concept of the robustness radius.

######  Definition 3.

For all nominal sequences $x\in\mathcal{X}_{k}^{*}$, let us consider the
robustness radius for the $k^{\prime}$th class label as

$\varepsilon_{k}=\min_{x\in\mathcal{X}_{k}^{*},\,j\neq
k}d\big{(}\,p(T;x),V_{j}\,\big{)},$ (15)

where $j=1,\dots,m$.

The robust radius quantifies the minimal distance from a nominal belief vector
$p(T;x)$ to the boundary of $V_{k}$, which enables us determine when the
classification with $\tilde{x}\in\mathcal{X}_{k}$ is robust, i.e., the result
is the same with the one generated by some $x\in\mathcal{X}_{k}^{*}$.

######  Theorem 2.

Suppose that the RNN model (2) is stable. The classification generated by a
perturbed sequence $\tilde{x}\in\mathcal{X}_{k}$ is robust for the $k$’th
class, if there exist some $x\in\mathcal{X}_{k}^{*}$ such that

$\|x-\tilde{x}\|_{\ell_{\infty}}<\varepsilon_{k}\,\eta^{-1},$ (16)

where $\|x-\tilde{x}\|_{\ell_{\infty}}$ and $\eta$ are defined in (8) and (1).

###### Proof.

In the view of the $k$’th class, the classification with $\tilde{x}\in X$ is
robust if

$d(\,p(T;\tilde{x}),V_{j}\,)=\inf_{p^{\prime}\in
V_{j}}d(p(T;\tilde{x}),p^{\prime})>0,$ (17)

for all $j=1,\dots,m$ and $j\neq k$. The above quantity obtains a lower bound
using the triangle inequality,

$\displaystyle\inf_{p^{\prime}\in V_{j}}d(p(T;\tilde{x}),p^{\prime})$
$\displaystyle\geq$ $\displaystyle\inf_{p^{\prime}\in V_{j}}$ $\displaystyle
d(p(T;x),p^{\prime})-d(p(T;x),p(T;\tilde{x})),$

for any $\tilde{x}\in\mathcal{X}_{k}$ and $x\in\mathcal{X}_{k}^{*}$. Then, the
inequality (17) will be satisfied if

$\displaystyle\inf_{p^{\prime}\in
V_{j}}d(\,p(T;x),p^{\prime}\,)-d(\,p(T;x),p(T;\tilde{x})\,)>0,$

which is equivalent to $\|p(T;x)-p(T;\tilde{x})\|<\varepsilon_{k}.$ Then, we
can conclude by applying Theorem 1 to the inequality above. ∎

The above theorem states that in the classification model with stable RNNs, a
perturbed data sequence $\tilde{x}$ will generate a correct classification
result if there exist some nominal sequences $x\in\mathcal{X}_{k}^{*}$ that
satisfy (16). In other words, the constant $\varepsilon_{k}\eta^{-1}$ measures
maximal allowed deviation for the input to perform robust classification on
the $k$’th class, i.e., for any perturbed sequence $\tilde{x}$ that satisfies
(16), the classification result is guaranteed to be correct. In Fig. 6, we
present the idea of this relation with a simplified example. We highlight that
the above result provides a sufficient condition, so $\tilde{x}$ may still
generate a correct classification if (16) is not satisfied.

Figure 6: This figure depicts nominal belief vectors in the $k$’th class
$p(T;x_{1}),...,p(T;x_{5})$, and the perturbed belief vector $p(T;\tilde{x})$.
The model will perform a robust classification if $p(T;\tilde{x})$ is within
the distance of $\varepsilon_{k}$ to some $p(T;x_{i})$, where
$x_{i}\in\mathcal{X}_{k}^{*}$.

##  VII Case Study and Simulations

We use the map and image classification [16] as examples to demonstrate and
validate our theoretical results. In the case study, it is assumed that an
aerial robot with a downward-facing camera aims to classify the underlying
image (or map). However, due to the limited sensing capability, the robot can
only observe a localized portion of the image at each time step, and it is
allowed to traverse the environment to collect partial images as a time-
indexed sequence.

###  VII-A Training for Classification

We consider the sampling routine for the robot is given and fixed for both
training and testing stage 222For each map (or image), we generate five unique
paths for a robot to traverse, and they are fixed through the training and
testing. However, robots are also capable of planning the path based on their
observation, see [17, 21]. Data sequences are generated with a fixed length
$T$. In the case study, a VGG-19 [27] model is adopted to process the raw
observations $O(t)\in\mathcal{O}$, such that $x(t)=VGG\big{(}O(t)\big{)}.$ The
robot also uses a stable LSTM cell [8], which is a special case of RNN, to
recurrently process the input. A stable (contractive) LSTM can be learned by
introducing the following constraints in the training stage 333In our
notation, $\|f\|_{\infty}=\sup_{t}\|f_{t\geq 0}\|_{\infty}$ and $f_{t}$ is the
output of the forget gate of the LSTM cell. We refer to [20] for the details.

$\max\left\\{\|W_{u}\|_{\infty},\|W_{o}\|_{\infty},4\|W_{z}\|_{\infty},\sqrt{\|W_{f}\|_{\infty}}\right\\}<1-\|f\|_{\infty}.$
(18)

When the above constraint is not satisfied during the training, each row of
the LSTM weight matrices will be divided by a constant scaling factor $\tau>1$
after each gradient step until (18) is satisfied.

In the training stage, we evaluate the classification reward with the data
sequence $x\in\mathcal{X}_{k}^{\prime}$ by a log-sum-exp (LSE) loss as
$r=-LSE\big{(}e_{k},p(T;x)\big{)},$ in which $k$ denotes the ground truth
label.

###  VII-B MNIST Dataset

Figure 7: This figure depicts the training progress with and without the
stability constraint (18) on the MNIST dataset.

In the first case study, we consider the robot is traveling over the image
from the MNIST dataset [14, 22]. The training and testing is performed in
PyTorch [23] with ADAM [10] and a learning rate $l_{r}=0.0001$. The testing
result is presented in Table VII-B, which is validated over five random seeds.

To establish a benchmark, we train an independent VGG-19 model with the entire
image (map) as the input, shown in the last column of Table VII-B. In the
second column, the observation size denotes the maximum possible coverage of
the entire image or map. It is shown that a single robot can classify by only
revealing a small portion of the environment. It should be emphasized that the
performance on the MNIST will reach $98\%$ by using multiple communicating
robots [22].

Dataset | Observation Size ($\%$) | Stability Constraint | Scaling Factor $\tau$ | Classification Accuracy ($\%$) | Benchmark (VGG-19)
---|---|---|---|---|---
Map Dataset with clouds | $\leq$7.8 | No | - | 77.62 | 99.43
Yes | 1.01 | 77.23
MNIST | $\leq$ 68.88 | No | - | 94.56 | 99.33
Yes | 1.05 | 94.34
Yes | 1.1 | 94.90
Yes | 2 | 95.05
Yes | 4 | 94.20
Yes | 8 | 94.04

TABLE I: The performance measured with different observation size, stability
condition and dataset.

###  VII-C Robustness Analysis on the MNIST Dataset

To reveal how a stable RNN model will affect the performance of the
classifier, we test both unstable and stable models with various scaling
factors $\tau$ in the MNIST dataset. As shown in both Fig. 7 and Table VII-B,
a stable model has a comparable performance and sometimes it even outperforms
the unstable model for certain values of $\tau$.

####  VII-C1 Constant deviation

To investigate the robustness of the classifier, we first obtain the nominal
sequences444A sequence is called nominal if its resulting classification
accuracy is $100\%$. from the training set $\mathcal{X}_{k}^{\prime}$ and
measure the performance of the model by adding a constant deviation
$\xi\in[0,10]$ to the nominal sequences to get ~x(t) = x(t) + ξ1_a, where
$\bm{1}_{a}$ is the vector of all ones. This implies that
$\|x-\tilde{x}\|_{\ell_{\infty}}=\xi$. The test is performed on the $5$’th
class of the MNIST dataset for the unstable and stable models with
$\tau=1.05$. The result is shown in Fig. 8, and it implies that, for the
stable model, if the deviation of the input $\|x-\tilde{x}\|_{\ell_{\infty}}$
is less than $\varepsilon_{5}/\eta=0.0886$, the classification with
$\tilde{x}$ is guaranteed to be robust, i.e., with a $100\%$ accuracy. This
agrees with our theoretical result. However, the unstable model starts to
generate wrong predictions before the quantity
$\|x-\tilde{x}\|_{\ell_{\infty}}$ reaches $0.0886$.

We highlight that for the stable model the accuracy remains $100\%$ even for
some $\|x-\tilde{x}\|_{\ell_{\infty}}>\varepsilon_{5}/\eta$. This is because
our theoretical result provides a sufficient condition, which is usually
conservative. As we observe from the simulations, there is still a chance to
get robust classification when (16) is not satisfied. Furthermore, it is
interesting to notice that the unstable model starts to outperform the stable
model when $\|x-\tilde{x}\|_{\ell_{\infty}}\geq 0.2$. The reason is that by
imposing the stability constraint (18) one only require the stable model to
confidently exhibit the robust classification with the input perturbation less
than $\varepsilon_{k}/\eta$, instead of concerning the performance with
$\|x-\tilde{x}\|_{\ell_{\infty}}\geq\varepsilon_{k}/\eta$. The stable model
provides the confidence of $100\%$ accuracy with all deviations
$\|x-\tilde{x}\|_{\ell_{\infty}}<0.0886$, which is not guaranteed for the
unstable models, as shown in the shaded area in Fig. 8. This difference is
crucial when the robot is performing high-precision tasks, in which any level
of mistake is not acceptable. On the other hand, the stable model that focuses
on improving the performance with the Gaussian input noise has been proposed
and validated in our previous work [1].

####  VII-C2 Variable deviation

It is also interesting to see how the model will perform with the variable
deviations instead of the constant ones. The differences from the testing
dataset to the training dataset can be naturally considered variable
deviations since the variation of handwritten digits from the training to the
testing set can not be modeled as the constant deviation. Hence, we compare
the performances of different models on both training and testing sets, see
Fig. 9 and Table VII-C2, in which we denote the training set as “unperturbed”
and the testing set as “perturbed.”

The results of variable deviation show that the performance loss from
unperturbed data to the perturbed data is significantly reduced by using a
stable model. In addition, the stable model sometimes outperforms the unstable
model if a proper scaling factor is selected, e.g., $\tau=1.1$ and $\tau=2$.
We also highlight that there exists an potential trade-off between the
performance loss and the overall performance: As shown in Table VII-C2, a
higher scaling factor $\tau$ usually implies a more stable RNN and less
performance loss, but also a weaker overall performance since the stability
constraint will potentially drive the model away from the optimal classifier.

Figure 8: This figure depicts the classification accuracy when the constant
perturbation has been added to the nominal sequences in the $5$’th class in
the MNIST dataset with $\varepsilon_{5}/\eta=0.0886$. From left to right are
magnified details. The accuracy eventually converges to $10\%$ since the total
number of labels is $m=10$, and a random guess has $1/10$ of the chance to get
the correct label.

Stability | $\tau$ | Unperturbed | Perturbed | Performance
---|---|---|---|---
Constraint | Sequence | Sequence | Loss
No | - | 95.88 | 94.56 | 1.3163
Yes | 1.05 | 94.68 | 94.34 | 0.3396
Yes | 1.1 | 95.24 | 94.90 | 0.3416
Yes | 2 | 95.50 | 95.05 | 0.4467
Yes | 4 | 94.40 | 94.20 | 0.1978
Yes | 8 | 94.14 | 94.04 | 0.1062

TABLE II: The performance loss under perturbation with various scaling factor
$\tau$ on the MNIST dataset.

Figure 9: The above figures depict performance with both unperturbed and
perturbed data for all models trained with various scaling factor $\tau$ on
the MNIST dataset. Each data point represents the average result from $20000$
sampled data sequences.

###  VII-D Campus Map Dataset

In the second case study, we test our model on the Campus Map dataset [17].
The training and testing are performed in the same platform with the first
case study with a learning rate $l_{r}=0.00001$. The testing result is
presented in Table VII-B and Fig. 10, which is validated over five random
seeds. It is shown that the stable model obtains a comparable performance as
the unstable model on the Campus map dataset and enjoys the robustness
guarantee when the input perturbation satisfies certain conditions. It should
also be emphasized that the performance of the map classification will reach
$97\%$ by using a team of communicating robots [16].

We also successfully performed a real-world experiment of map classification
with aerial robots. Some snapshots taken from the experiments are shown in
Fig. 11, and a full experiment video can be found at
https://youtu.be/nsnPFAvJLoY.

Figure 10: This figure depicts the training progress with and without the
stability constraint (18) on the Campus map dataset. Figure 11: The above
figures are snapshots taken while the aerial robot is flying over the map of
the Lehigh University campus. The robot’s location and observation are denoted
with the red and the yellow box. The bar plot shows the real-time belief
vector. Eventually, the robot correctly classifies the map as the Lehigh
University (gray bar). For more detail, please see the full experiment video
at https://youtu.be/nsnPFAvJLoY.

##  VIII Conclusion

We present a framework to analyze the robustness properties of stable RNNs
with sequential inputs for classification purposes. It is shown that every
trained RNN exhibits robust classification with respect to some bounded
perturbations. We quantify robustness bounds in terms of trainable weight
matrices. Our results are significant as they reveal interplay among various
design (trainable) parameters. Our extensive simulations and one real-world
experiment support and validate the usefulness of our theoretical findings.

## References

  * [1] Arash Amini, Guangyi Liu and Nader Motee  “Robust Learning of Recurrent Neural Networks in Presence of Exogenous Noise”  In _2021 60th IEEE Conference on Decision and Control (CDC)_ , 2021, pp. 783–788
  * [2] M. Bock, Jürgen Böhner, Olaf Conrad, Rüdiger Köthe and Andre Ringeler  “XV. Methods for creating Functional Soil Databases and applying Digital Soil Mapping with SAGA GIS”  In _JRC Scientific and technical Reports, Office for Official Publications of the European Communities, Luxemburg_ , 2007
  * [3] R.. Boie and I.. Cox  “An analysis of camera noise”  In _IEEE Transactions on Pattern Analysis & Machine Intelligence_ 14.06 IEEE Computer Society, 1992, pp. 671–674
  * [4] A. Breitenmoser, M. Schwager, Jean-Claude Metzger, R. Siegwart and D. Rus  “Voronoi coverage of non-convex environments with a group of networked robots”  In _2010 IEEE international conference on robotics and automation_ , 2010  IEEE
  * [5] Nicholas Carlini, Anish Athalye, Nicolas Papernot, Wieland Brendel, Jonas Rauber, Dimitris Tsipras, Ian Goodfellow, Aleksander Madry and Alexey Kurakin  “On evaluating adversarial robustness”  In _arXiv preprint arXiv:1902.06705_ , 2019
  * [6] F. De la Torre and M.. Black  “Robust principal component analysis for computer vision”  In _Proceedings Eighth IEEE International Conference on Computer Vision. ICCV 2001_ 1, 2001, pp. 362–369  IEEE
  * [7] Noel E Du Toit and Joel W Burdick  “Robot motion planning in dynamic, uncertain environments”  In _IEEE Transactions on Robotics_ 28.1 IEEE, 2011, pp. 101–115
  * [8] S. Hochreiter and J. Schmidhuber  “Long Short-term Memory”  In _Neural computation_ 9, 1997, pp. 1735–80
  * [9] Thomas S Huang and RY Tsai  “Image sequence analysis: Motion estimation”  In _Image sequence analysis_ Springer, 1981
  * [10] D.. Kingma and J. Ba  “Adam: A Method for Stochastic Optimization”, 2014  arXiv:1412.6980 [cs.LG]
  * [11] R. Klein  “Abstract Voronoi diagrams and their applications”  In _Workshop on Computational Geometry_ , 1988  Springer
  * [12] Ching-Yun Ko, Zhaoyang Lyu, Lily Weng, Luca Daniel, Ngai Wong and Dahua Lin  “POPQORN: Quantifying robustness of recurrent neural networks”  In _International Conference on Machine Learning_ , 2019, pp. 3468–3477  PMLR
  * [13] D. Kohli  “Machine Learning: Is the softmax function Lipschitz with Lipschitz constant 1?” URL:https://math.stackexchange.com/q/2021011, Mathematics Stack Exchange
  * [14] Y. LeCun, Léon Bottou, Yoshua Bengio and Patrick Haffner  “Gradient-based learning applied to document recognition”  In _Proceedings of the IEEE_ 86.11 Ieee, 1998
  * [15] C. Li and S. Zhang  “Stationary probability vectors of higher-order Markov chains”  In _Linear Algebra and Its Applications_ 473 Elsevier, 2015, pp. 114–125
  * [16] G. Liu, A. Amini, M. Takáč, H. Muñoz-Avila and N. Motee  “Distributed Map Classification using Local Observations”, 2020  arXiv:2012.10480 [cs.RO]
  * [17] Guangyi Liu, Arash Amini, Martin Takáč and Nader Motee  “Classification-Aware Path Planning of Network of Robots”  In _International Symposium Distributed Autonomous Robotic Systems_ , 2021, pp. 294–305  Springer
  * [18] Daniel Lowd and Christopher Meek  “Adversarial learning”  In _Proceedings of the eleventh ACM SIGKDD international conference on Knowledge discovery in data mining_ , 2005, pp. 641–647
  * [19] Alexander Matyasko and Lap-Pui Chau  “Improved network robustness with adversary critic”  In _Advances in Neural Information Processing Systems_ 31, 2018
  * [20] J. Miller and M. Hardt  “Stable Recurrent Models”  In _International Conference on Learning Representations_ , 2018
  * [21] H.. Mousavi, Guangyi Liu, Weihang Yuan, Martin Takáč, Héctor Muñoz-Avila and Nader Motee  “A Layered Architecture for Active Perception: Image Classification using Deep Reinforcement Learning”, 2019  arXiv:1909.09705 [cs.LG]
  * [22] H.. Mousavi, Mohammadreza Nazari, Martin Takáč and Nader Motee  “Multi-Agent Image Classification via Reinforcement Learning”  In _2019 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS)_ , 2019
  * [23] A. Paszke, Sam Gross, Soumith Chintala, Gregory Chanan, Edward Yang, Zachary DeVito, Zeming Lin, Alban Desmaison, Luca Antiga and Adam Lerer  “Automatic differentiation in pytorch”, 2017
  * [24] Visvanathan Ramesh and Robert M Haralick  “Random perturbation models and performance characterization in computer vision”  In _Proceedings 1992 IEEE Computer Society Conference on Computer Vision and Pattern Recognition_ , 1992, pp. 521–522  IEEE Computer Society
  * [25] Alexander Robey, Haimin Hu, Lars Lindemann, Hanwen Zhang, Dimos V Dimarogonas, Stephen Tu and Nikolai Matni  “Learning control barrier functions from expert demonstrations”  In _2020 59th IEEE Conference on Decision and Control (CDC)_ , 2020  IEEE
  * [26] W. Rudin  “Real and Complex Analysis”  McGraw-Hill, 1974
  * [27] K. Simonyan and A. Zisserman  “Very deep convolutional networks for large-scale image recognition”  In _3rd International Conference on Learning Representations_ , 2015
  * [28] Justin Solomon  “Numerical algorithms: methods for computer vision, machine learning, and graphics”  CRC press, 2015
  * [29] Niko Sünderhauf, Oliver Brock, Walter Scheirer, Raia Hadsell, Dieter Fox, Jürgen Leitner, Ben Upcroft, Pieter Abbeel, Wolfram Burgard and Michael Milford  “The limits and potentials of deep learning for robotics”  In _The International Journal of Robotics Research_ 37.4-5 SAGE Publications Sage UK: London, England, 2018, pp. 405–420
  * [30] Zhongwei Tang, Rafael Grompone Gioi, Pascal Monasse and Jean-Michel Morel  “A precision analysis of camera distortion models”  In _IEEE Transactions on Image Processing_ 26.6 IEEE, 2017
  * [31] Kaidi Xu, Sijia Liu, Pu Zhao, Pin-Yu Chen, Huan Zhang, Quanfu Fan, Deniz Erdogmus, Yanzhi Wang and Xue Lin  “Structured adversarial attack: Towards general implementation and better interpretability”  In _arXiv preprint arXiv:1808.01664_ , 2018