# Temporal Logic Resilience for Dynamical Systems⋆

Adnane Saoud1 ,  Pushpak Jagtap2 and  Sadegh Soudjani3 1College of Computing,
University Mohammed VI Polytechnic, Benguerir, Morocco<EMAIL_ADDRESS>2Robert Bosch Center for Cyber-Physical Systems, Indian Institute of Science,
Bangalore, India<EMAIL_ADDRESS>3Max Planck Institute for Software
Systems, Kaiserslautern, D-67663, Germany<EMAIL_ADDRESS>
###### Abstract.

We consider the notion of resilience for cyber-physical systems, that is, the
ability of the system to withstand adverse events while maintaining acceptable
functionality. We use finite temporal logic to express the requirements on the
acceptable functionality and define the resilience metric as the maximum
disturbance under which the system satisfies the temporal requirements. We fix
a parameterized template for the set of disturbances and form a robust
optimization problem under the system dynamics and the temporal specifications
to find the maximum value of the parameter. Additionally, we introduce two
novel classes of specifications: closed and convex finite temporal logics
specifications, offering a comprehensive analysis of the resilience metric
within these specific frameworks. From a computational standpoint, we present
an exact solution for linear systems and exact-time reachability and finite-
horizon safety, complemented by an approximate solution for finite-horizon
reachability. Extending our findings to nonlinear systems, we leverage linear
approximations and SMT-based approaches to offer viable computational
methodologies. The theoretical results are demonstrated on the temperature
regulation of buildings, adaptive cruise control and DC motors.

⋆ This work was supported in part by by the ARTPARK, the Google Research
Grant, the SERB Start-up Research Grant, and the CSR Grants by Siemens and
Nokia. This project is also supported by the following grants: EPSRC
EP/V043676/1, EIC 101070802, and ERC 101089047.

##  1\. Introduction

Resilience has been studied by many research communities and it is broadly
defined as _the ability of a system to withstand adverse events while
maintaining an acceptable functionality._ For critical infrastructures,
resilience is the main factor determining their reliability and is improved by
continuously enhancing the prevention and absorption of disruptive events, and
the recovery and adaptation for such events [RSHL19]. For IT systems,
resilience is considered mainly against adverse cyber events, which are the
cyber attacks that negatively impact the availability, integrity, or
confidentiality of the system [Pil19]. With climate change increasing the
extreme flood events, resilience metrics that consider the dynamical changes
of the system have also received attention in the water research community to
define and assess the resilience of water resource recovery facilities
[HWOY22].

In this paper, we provide a notion of resilience for cyber-physical systems
that integrates the time-evolution of the system with temporal logic to
provide a quantitative measure on how the system can cope with disturbances.
The temporal logic is used to formally encode the safety and other compliance
requirements on the operation of the system and also express the expected
behavior of the system to disturbances. We define resilience as the largest
disturbance within a given (parameterized) set that can be applied to the
system in its time evolution while still satisfying the temporal logic
specification. This can also be interpreted as the minimum disturbance that
needs to be applied to the system to falsify the specification.

We focus on a discrete-time dynamical model of the system and express
resilience requirements as linear temporal logic specifications over finite
traces (LTLF) [ZTL+17]. We provide different characterizations of the
resilience metrics for general LTLF specifications. Moreover, we introduce two
novel classes of LTLF specifications: closed and convex finite LTLF
specifications. We then characterize the fragments of LTLF specifications that
are closed and convex, and provide an analysis of the structural properties of
resilience with respect to closed and convex specifications. From a
computation perspective, we show how the optimization for computing resilience
can be solved exactly for linear systems with respect to specifications,
including the exact-time reachability and finite-horizon safety, and
approximately with respect to finite-horizon reachability. Moreover, we
provide under-approximations of the resilience metric for nonlinear systems
using linear programs and SMT-based approaches. We then provide numerical
examples showing the merits of the proposed results.

Our definition of resilience is substantially different from _robustness_ for
temporal specifications, which is defined as follows. Robustness of a system
$\Sigma$ with respect to a temporal specification $\phi$ is the largest value
$\varepsilon$ such that we still satisfy $\phi$ if we expand the solutions of
$\Sigma$ with a uniform $\varepsilon$ over time and over trajectories [DM10,
FP09]. This definition does not take into account the dynamics of the system
$\Sigma$ and directly applies the expansion to the solution of $\Sigma$. In
reality, disturbances and extreme events affect the time evolution of the
system, and this needs to be integrated with any definition of resilience for
dynamical systems.

Related work. The literature on defining quantitative semantics for different
classes of temporal logic is relatively mature. These quantitative semantics
study how well the system trajectories satisfy a given specification. The
techniques include using discounting modalities that give less importance to
distant events [ABK14] and averaging modalities [BMM14] where the semantics of
standard modalities are extended using min, max, and a long-run average
operator. The paper [DM10] considers real-valued signals and presents variants
of robustness measures that indicate how far a given signal stands, in space
and time, from satisfying or violating a property and studies their
sensitivity to the parameters of the system. The paper [FP09] considers the
robust interpretation of Metric Temporal Logic to connect robust satisfaction
of properties on discrete-time signals to their continuous-time counterparts.

The authors in [DFM13] present an efficient algorithm for computing the
robustness degree in which a piecewise-continuous signal satisfies or violates
a Signal Temporal Logic (STL) formula. Application of robustness metrics in
specification-based monitoring of cyber-physical systems (CPS) is provided in
[BDD+18] with a survey of theory and tools. The robustness metric is also used
for temporal logic falsification of CPS [ATMMR18]. STL is also used in
[CLSP22, CSPL23] to study two important resilience properties of CPS, which
are recoverability and durability. Having tolerance thresholds in both time
and value is considered in [GM20] to characterize approximate notions of
hybrid conformance. The notion of robustness with respect to a temporal
specification is also defined and used in stochastic systems. The paper
[IHL23] defines robustness for continuous-time linear stochastic systems with
respect to STL specifications. This notion computes the probability of
satisfying atomic prepositions in the STL formula and combines the related
probabilities using logical connectives based on the structure of the STL
formula to get an interval for robustness. The robustness of STL
specifications is used in [FMPS18] to perform controller synthesis for
stochastic systems using a model predictive control approach. While we do not
consider stochastic disturbances in this work, any available statistics on the
likelihood of disturbances can be mapped to statistical properties of the
resilience metric, e.g., by constructing a distribution for the resilience
value when the distribution of the possible disturbances is known.

All the works mentioned above study the robust satisfaction of properties for
a given set of disturbances. In contrast, we are looking at characterizing
resilience to compute the largest disturbance within a given parameterized
set. The work closest in spirit to our approach is the paper [SDSC17] that is
limited to linear systems and safety specifications and studies properties of
polytopic invariant sets when the size of the disturbance set changes with a
scaling factor. The current paper also extends the preliminary results
presented in [SJS23]. In the current version, we added numerical examples and
a new case study to illustrate the proposed concepts, we provided a new
approach to compute the resilience metric for finite-horizon reachability
using scenario optimization [ESL14] and a new approach to compute the
resilience metric for nonlinear systems using an SMT-based approach [BT18].
Finally, we presented a characterization of the fragments of LTLF
specifications that are closed and convex, and we included all the proofs.

In summary, the main contributions of this paper are as follows. We provide a
quantitative resilience metric by integrating the underlying dynamics with
temporal logic specifications. We show how the related optimizations can be
solved for linear systems and various types of specifications. We show that
resilience with respect to convex or closed specifications enjoys some _nice_
properties. Finally, we provide under-approximations of the resilience metric
for nonlinear systems.

The rest of this paper is organized as follows. We provide the preliminaries
on the model class and temporal logic specifications in Sec. 2. We define the
notion of resilience and study its structural properties in Sec. 3. We provide
computations for linear systems in Sec. 4 and discuss the computation of lower
bounds for nonlinear systems in Sec. 5. The results are demonstrated on
examples from adaptive cruise control, temperature regulation of buildings,
and DC motors in Sec. 6 with concluding remarks in Sec. 7.

##  2\. Preliminaries

Notation: The symbols $\mathbb{R}$, $\mathbb{R}_{\geq 0}$, $\mathbb{N}$, and
$\mathbb{N}_{\geq n}$ denote the set of real, nonnegative real numbers,
nonnegative integers, and integers greater than or equal to $n\in\mathbb{N}$,
respectively. We use $\mathbb{R}^{n\times m}$ to denote the space of real
matrices with $n$ rows and $m$ columns. For a matrix $A\in\mathbb{R}^{n\times
m}$, $A^{T}$ represent the transpose of $A$. For a vector
$x\in\mathbb{R}^{n}$, we use $\|x\|$ and $\|x\|_{\infty}$ to denote the
Euclidean and infinity norm, respectively. We use $\mathbb{I}$ to denote the
identity matrix. For a set of $p$ points $C=\\{c_{1},c_{2},\ldots,c_{p}\\}$,
$c_{i}\in\mathbb{R}^{n}$, the convex hull of $C$ is represented by
$\textrm{conv}(c_{1},c_{2},\ldots,c_{p}):=\\{\alpha_{1}c_{1}+\alpha_{2}c_{2}+\ldots+\alpha_{p}c_{p}\mid
c_{i}\in C,\alpha_{i}\geq
0,i\in\\{1,2,\ldots,p\\},\sum_{i=1}^{p}\alpha_{i}=1\\}$. An interval in
$\mathbb{R}^{n}$ is a set denoted by
$X=[\underline{x}_{1},\overline{x}_{1}]\times[\underline{x}_{2},\overline{x}_{2}]\times\ldots\times[\underline{x}_{n},\overline{x}_{n}]$
and defined as $X=\left\\{x\in\mathbb{R}^{n}\mid\underline{x}_{i}\leq
x_{i}\leq\overline{x}_{i},~{}i\in\\{1,2,\ldots,n\\}\right\\}$, where
$x_{i}\in\mathbb{R}$ represents the ith component of the vector
$x\in\mathbb{R}^{n}$. In particular, when $\underline{x}_{i}=\underline{x}$
and $\overline{x}_{i}=\overline{x}$ for all $i\in\\{1,2,\ldots,n\\}$, then the
interval $X$ can be written in a compact form as:
$X=[\underline{x},\overline{x}]^{n}$. Given $x\in\mathbb{R}^{n}$ and
$\varepsilon\geq 0$,
$\Omega_{\varepsilon}(x)=\\{z\in\mathbb{R}^{n}\mid\|z-x\|_{\infty}\leq\varepsilon$}.

###  2.1. Discrete-Time Dynamical Systems

A discrete-time system is a tuple $\Sigma=(X,D,f)$, where
$X\subset\mathbb{R}^{n}$ is the state space and $D\subset\mathbb{R}^{n}$ is
the disturbance space which is assumed to be a compact set containing the
origin. The evolution of the state of $\Sigma$ is given by

(1) $x(k+1)=f(x(k))+d(k),~{}~{}k\in\mathbb{N},$

where $d(k)\in D$ represents the additive disturbance. The trajectory of
system $\Sigma$ of length $N$ is represented by
$w_{x}=(x_{0},x_{1},\ldots,x_{N-1})\in X^{N}$, where $x_{k}$ represents the
value of trajectory starting from a state $x(0)=x_{0}\in X$ at $k^{\text{th}}$
instance (i.e., $x(k)$).

###  2.2. LTLF Specifications

Linear temporal logic (LTL) provides a high-level language for describing the
desired behavior of a dynamical system. Formulas in this logic are constructed
inductively using a set of atomic propositions and combining them via Boolean
operators [BK08]. Consider a finite set of atomic propositions $\mathsf{AP}$
that defines the alphabet $\Sigma_{a}:=2^{\mathsf{AP}}$. Each letter of this
alphabet evaluates a subset of the atomic propositions as true. In this work,
we consider LTL specifications over finite words, referred to as LTLF, where
the letters form finite words defined as
$\boldsymbol{\omega}=(\boldsymbol{\omega}_{0},\boldsymbol{\omega}_{1},\boldsymbol{\omega}_{2},\ldots,\boldsymbol{\omega}_{N-1})\in\Sigma_{a}^{N}$
for some $N\in\mathbb{N}$, with $\boldsymbol{\omega}_{i}$ representing the
letter of the word at $i$th instance. These words are connected to
trajectories of the system $\Sigma$ via a measurable labeling function
$\textsf{L}:X\rightarrow\Sigma_{a}$ that assigns letters
$\alpha=\textsf{L}(x)$ to state $x\in X$. That is, any finite trajectory
$w_{x}=(x_{0},x_{1},\ldots,x_{N-1})$ is mapped to the set of finite words
$\Sigma_{a}^{N}$, as
$\boldsymbol{\omega}=\textsf{L}(w_{x}):=(\textsf{L}(x_{0}),\textsf{L}(x_{1}),\textsf{L}(x_{2}),\ldots\textsf{L}(x_{N-1}))$.

######  Definition 2.1.

An LTLF formula over a set of atomic propositions $\mathsf{AP}$ is constructed
inductively as

$\psi::=\textsf{true}\,|\,p\,|\,\neg\psi\,|\,\psi_{1}\wedge\psi_{2}\,|\,\psi_{1}\vee\psi_{2}\,|\,\mathord{\bigcirc}\psi\,|\,\psi_{1}\mathbin{\sf
U}\psi_{2}\,|\,\square\psi\,|\,\lozenge\psi,$

with $p\in\mathsf{AP}$, and $\psi_{1},\psi_{2},\psi$ being LTLF formulas.

Given a finite word $\boldsymbol{\omega}$ of length $N$ and an LTLF formula
$\psi$, we inductively define when an LTLF formula is true at the
$i^{\text{th}}$ step $(i<N)$ and denoted by
$\boldsymbol{\omega}_{i}\models\psi$, as follows:

  * •

$\boldsymbol{\omega}_{i}\vDash\textsf{true}$ always hold and
$\boldsymbol{\omega}_{i}\vDash\textsf{false}$ does not hold.

  * •

An atomic proposition, $\boldsymbol{\omega}_{i}\vDash p$ for $p\in\mathsf{AP}$
holds if $p\in\boldsymbol{\omega}_{i}$.

  * •

A negation, $\boldsymbol{\omega}_{i}\vDash\neg p$, holds if
$\boldsymbol{\omega}_{i}\nvDash p$.

  * •

A logical conjunction, $\boldsymbol{\omega}_{i}\vDash\psi_{1}\wedge\psi_{2}$,
holds if $\boldsymbol{\omega}_{i}\vDash\psi_{1}$ and
$\boldsymbol{\omega}_{i}\vDash\psi_{2}$.

  * •

A logical disjunction, $\boldsymbol{\omega}_{i}\vDash\psi_{1}\vee\psi_{2}$,
holds if $\boldsymbol{\omega}_{i}\vDash\psi_{1}$ or
$\boldsymbol{\omega}_{i}\vDash\psi_{2}$.

  * •

A temporal next operator,
$\boldsymbol{\omega}_{i}\vDash\mathord{\bigcirc}\psi$, holds if
$\boldsymbol{\omega}_{i+1}\vDash\psi$. Similarly, for $0\leq j<N-i$,
$\boldsymbol{\omega}_{i}\vDash\mathord{\bigcirc}^{j}\psi$, holds if
$\boldsymbol{\omega}_{i+j}\vDash\psi$.

  * •

A temporal until operator, $\boldsymbol{\omega}_{i}\vDash\psi_{1}\mathbin{\sf
U}\psi_{2}$, holds if for some $m$ such that $i\leq m<N$, we have
$\boldsymbol{\omega}_{m}\vDash\psi_{2}$ and for all $i\leq k<m$, we have
$\boldsymbol{\omega}_{k}\vDash\psi_{1}$.

  * •

A temporal always operator, $\boldsymbol{\omega}_{i}\vDash\square\psi$, holds
if for all $m$ such that $i\leq m<N$, we have
$\boldsymbol{\omega}_{m}\vDash\psi$. Similarly, for $0\leq j<N-i$,
$\boldsymbol{\omega}_{i}\vDash\square^{j}\psi$ holds if for all $i\leq m\leq
i+j$, we have $\boldsymbol{\omega}_{m}\vDash\psi$.

  * •

A temporal eventually operator, $\boldsymbol{\omega}_{i}\vDash\lozenge\psi$,
holds if for some $m$ such that $i\leq m<N$, we have
$\boldsymbol{\omega}_{m}\vDash\psi$. Similarly, for $0\leq j<N-i$,
$\boldsymbol{\omega}_{i}\vDash\lozenge^{j}\psi$, holds if for some $m$ such
that $i\leq m<i+j$, we have $\boldsymbol{\omega}_{m}\vDash\psi$.

An LTLF formula $\psi$ is true for $\boldsymbol{\omega}$, denoted by
$\boldsymbol{\omega}\vDash\psi$, if and only if
$\boldsymbol{\omega}_{0}\vDash\psi$. For a trajectory
$w_{x}=(x_{0},x_{1},\ldots,x_{N-1})\in X^{N}$, we say that $w_{x}\vDash\psi$
if for
$\boldsymbol{\omega}=\textsf{L}(w_{x}):=(\textsf{L}(x_{0}),\textsf{L}(x_{1}),\ldots,\textsf{L}(x_{N-1}))$,
we have $\boldsymbol{\omega}\vDash\psi$. Similarly, for a set of trajectories
$\mathcal{X}\subseteq X^{N}$, we say that $\mathcal{X}\vDash\psi$, if
$w_{x}\vDash\psi$ for all $w_{x}\in\mathcal{X}$.

######  Remark 1.

Our notion of resilience is general, but for computational purposes, we
restrict ourselves to the following specifications over words of length $N$:
$\square p$, $\mathord{\bigcirc}p$, $\lozenge p$, with $p\in\mathsf{AP}$, and
conjunctions over them. Note that all LTLF can be represented using
Deterministic Finite Automata (DFA) [ZPV19] and effectively, one can represent
them using sequences of reach and avoid specifications (i.e., $\psi=\lozenge
p\wedge\square\neg q$, where $p,q\in\Sigma_{a}$) [WKK22, Section III.b]. Thus,
one can easily use the results provided in the paper for any LTLF
specification using the properties provided in Proposition 3.1.

##  3\. Characterizations of Resilience

The goal of the paper is to provide characterizations and algorithmic
procedures for computing resilience with respect to different classes of
systems and specifications.

###  3.1. Resilience for LTLF Specifications

Consider the system $\Sigma$ in (1), with a set of disturbances given by a
ball with respect to infinity norm centered at zero:
$D:=\Omega_{\varepsilon}(0)$. We denote by $\xi(x,\varepsilon)$ the set of
trajectories starting from some $x\in X$ with such a disturbance set:

(2)
$\xi(x,\varepsilon):=\left\\{(x_{0},x_{1},x_{2},\ldots)\,|\,x_{0}=x,\,x_{k+1}\in
f(x_{k})+D\right\\}.$

Note that $\xi(x,0)$ contains only the disturbance-free trajectory of the
system (nominal trajectory). Similarly, given a set $A\subseteq X$, we define
the set of trajectories $\xi(A,\varepsilon)=\bigcup\limits_{x\in
A}\xi(x,\varepsilon).$

######  Definition 3.1.

Consider the dynamical system $\Sigma$ in (1), an LTLF specification $\psi$
and a point $x\in X$. We define the resilience of the system $\Sigma$ with
respect to the initial condition $x$ and the specification $\psi$ as a
function $g_{\psi}:X\rightarrow\mathbb{R}_{\geq 0}\cup\\{+\infty\\}$ with:

(3) $g_{\psi}(x)=\begin{cases}\sup\left\\{\varepsilon\geq
0\,|\,\xi(x,\varepsilon)\vDash\psi\right\\},&\text{ if }\xi(x,0)\vDash\psi\\\
0,&\text{ if }\xi(x,0)\nvDash\psi,\end{cases}$

where $\xi(x,\varepsilon)\vDash\psi$ indicates that all trajectories in
$\xi(x,\varepsilon)$ satisfy the specification.

Similarly, for a set $A\subseteq X$, the resilience metric over the set $A$
and specification $\psi$ can be defined as follows:

(4) $g_{\psi}(x)=\begin{cases}\sup\left\\{\varepsilon\geq
0\,|\,\xi(x,A)\vDash\psi\right\\},&\text{ if }\xi(x,A)\vDash\psi\\\ 0,&\text{
if }\xi(x,A)\nvDash\psi,\end{cases}$

######  Remark 2.

The above definition assigns a non-negative value to resilience; hence we call
it a _resilience metric_. We can extend this definition to take both positive
and negative values depending on the satisfaction of the specification by the
nominal trajectory $\xi(x,0)$. For a given $\psi$, let us define the function
$h_{\psi}:X\rightarrow\mathbb{R}\cup\\{\pm\infty\\}$ with

(5) $h_{\psi}(x):=\begin{cases}\sup\left\\{\varepsilon\geq
0\,|\,\xi(x,\varepsilon)\vDash\psi\right\\},&\text{ if }\xi(x,0)\vDash\psi\\\
-\sup\left\\{\varepsilon\geq
0\,|\,\xi(x,\varepsilon)\vDash\neg\psi\right\\},&\text{ if
}\xi(x,0)\not\vDash\psi.\end{cases}$

Note that

(6) $h_{\psi}(x)=\begin{cases}g_{\psi}(x),&\text{ if }\xi(x,0)\vDash\psi\\\
-g_{\neg\psi}(x),&\text{ if }\xi(x,0)\not\vDash\psi.\end{cases}$

This definition is useful for guiding the system towards satisfying the
specification and then maximizing the resilience by appropriately quantifying
disturbances and decision variables. In this paper, we study the resilience of
autonomous systems (i.e., systems without control inputs or decision
variables). Therefore, we use Definition 3.1 in the rest of the paper.

###  3.2. Structural Properties of Resilience

In this subsection, we prove the structural properties of the resilience
metric in Definition. 3.1 utilizing the inductive definition of temporal
specifications.

We first have the following set of properties of the resilience.

######  Proposition 3.1.

Consider the dynamical system $\Sigma$ in (1), an LTLF specification $\psi$
and a point $x\in{X}$. The following properties hold:

  * (i)

When $\psi$ is the true specification, $g_{\psi}(x)=+\infty$, for all $x\in
X$.

  * (ii)

When $\psi$ is the false specification, $g_{\psi}(x)=0$, for all $x\in X$.

  * (iii)

For any specification $\psi=\psi_{1}\wedge\psi_{2}$, we have that
$g_{\psi}(x)=\min\\{g_{\psi_{1}}(x),g_{\psi_{2}}(x)\\}$, for all $x\in X$.

  * (iv)

For any specification $\psi=\psi_{1}\vee\psi_{2}$, we have that
$g_{\psi}(x)\geq\max\\{g_{\psi_{1}}(x),g_{\psi_{2}}(x)\\}$, for all $x\in X$.

  * (v)

For any specification $\psi=\neg\phi$, we have that
$g_{\psi}(x)g_{\phi}(x)=0$, for all $x\in X$.

  * (vi)

For $A\subset X$, $g_{\psi}(X)=\inf_{x\in X}g_{\psi}(x)$.

###### Proof.

Properties (i), (ii), and (vi) are the direct consequence of Definition 3.1.

Proof for property (iii): For $x\in X$, consider sets $A=\\{\epsilon\geq
0\mid\xi(x,\epsilon)\vDash\psi_{1}\\}$ and $B=\\{\epsilon\geq
0\mid\xi(x,\epsilon)\vDash\psi_{2}\\}$. Then we can write $A\cap
B=\\{\epsilon\geq 0\mid\xi(x,\epsilon)\vDash\psi_{1}\text{ and
}\xi(x,\epsilon)\vDash\psi_{2}\\}$. By the fact that sets $A$ and $B$ are
bounded from below, we have $\sup\\{A\cap B\\}=\min\\{\sup A,\sup B\\}$. This
implies
$g_{\psi_{1}\wedge\psi_{2}}(x)=\min\\{g_{\psi_{1}}(x),g_{\psi_{2}}(x)\\}$.

Proof for property (iv): For $x\in X$, consider sets $A=\\{\epsilon\geq
0\mid\xi(x,\epsilon)\vDash\psi_{1}\\}$ and $B=\\{\epsilon\geq
0\mid\xi(x,\epsilon)\vDash\psi_{2}\\}$. Then we can write $A\cup
B=\\{\epsilon\geq 0\mid\xi(x,\epsilon)\vDash\psi_{1}\text{ or
}\xi(x,\epsilon)\vDash\psi_{1}\\}$. By the fact that $\sup\\{A\cup
B\\}\geq\max\\{\sup A,\sup B\\}$, one has
$g_{\psi_{1}\vee\psi_{2}}(x)\geq\max\\{g_{\psi_{1}}(x),g_{\psi_{2}}(x)\\}$.

Proof for property (v): We consider two cases:
case (i): for $x\in X$, consider $\xi(x,0)\vDash\phi$, then
$g_{\psi=\neg\phi}(x)=0$. Thus $g_{\psi}(x)g_{\phi}(x)=0$.
case (ii): for $x\in X$, consider $\xi(x,0)\vDash\psi=\neg\phi$, then
$g_{\phi}(x)=0$. Thus $g_{\psi}(x)g_{\phi}(x)=0$.
This concludes the proof. ∎

We now provide an example showing that equality in the property (iv) does not
hold in general.

######  Example 3.1.

Consider the linear system $\Sigma$ with

$\displaystyle x_{1}(k+1)$ $\displaystyle=x_{1}(k)+3x_{2}(k)+d_{1}(k)$
$\displaystyle x_{2}(k+1)$ $\displaystyle=3x_{1}(k)-2x_{2}(k)+d_{2}(k)$

Consider the sets $S_{1}=[2,3]\times[1,2]$, $S_{2}=[2,3]\times[2,3]$ and
$S=S_{1}\cup S_{2}$. Consider the specifications
$\psi_{1}=\mathord{\bigcirc}S_{1}$, $\psi_{2}=\mathord{\bigcirc}S_{2}$ and
$\psi=\psi_{1}\vee\psi_{2}=\mathord{\bigcirc}(S_{1}\cup S_{2})$. For the state
$x(0)=(1,0.5)$, since $x(1)=(2.5,2)$ is both on the boundaries of the sets
$S_{1}$ and $S_{2}$, one can easily check that
$g_{\psi}(x)=0.5>\max(g_{\psi_{1}}(x),g_{\psi_{2}}(x))=0$.

###  3.3. Resilience Properties for Closed Specifications

Now, we provide sufficient conditions on the dynamics of the system $\Sigma$
and the specification $\psi$ allowing us to replace the $\sup$ operator with
the $\max$ operator in the definition of $g_{\psi}$ in (3), which makes the
computation of the resilience metric tractable. To do this, we introduce the
class of closed specifications defined below. Given a sequence of trajectories
$w_{x,i}$, $i\in\mathbb{N}$, with $w_{x,i}=x_{0,i},x_{1,i},x_{2,i}\ldots$, the
limit of $w_{x,i}$ is defined by
$w_{x}=\lim_{i\rightarrow\infty}w_{x,i}=x_{0},x_{1},x_{2}\ldots$, where for
all $j\in\mathbb{N}$, $x_{j}=\lim_{i\rightarrow\infty}x_{j,i}$. The sequence
$w_{x,i}$ is called converging whenever $w_{x}$ exists. Hence, the limit of a
sequence of trajectories can be seen as an element-wise limit of its
components.

######  Definition 3.2 (Closed specification).

Consider a set $X\subseteq\mathbb{R}^{n}$, an LTLF formula $\psi$ is said to
be closed if the following holds: for any converging sequence of trajectories
$w_{x,i}$, $i\in\mathbb{N}$, if for all $i\in\mathbb{N}$, $w_{x,i}\vDash\psi$,
then $w_{x}=\lim_{i\rightarrow\infty}w_{x,i}\vDash\psi$.

Intuitively, the closedness property states that a specification is preserved
when going from a sequence of trajectories to its element-wise limit. In
Appendix 7.2, we provide different characterizations of the fragment of LTLF
specifications that are closed.

Before stating the main result, we have the following auxiliary lemma.

######  Lemma 3.2.

Consider the discrete-time system $\Sigma$ in (1), an $\varepsilon\geq 0$ and
an LTLF formula $\psi$. If $\xi(x,\varepsilon)\vDash\psi$ then
$\xi(x,\varepsilon^{\prime})\vDash\psi$ for any
$\varepsilon^{\prime}\leq\varepsilon$.

###### Proof.

The proof follows directly from the fact that
$\xi(x,\varepsilon^{\prime})\subseteq\xi(x,\varepsilon)$ for
$\varepsilon^{\prime}\leq\varepsilon$. ∎

######  Proposition 3.3.

Consider the discrete-time system $\Sigma$ in (1) defined on a metric space
$X$. Consider $x\in X$ and an LTLF specification $\psi$. If the map
$f:\mathbb{R}^{n}\rightarrow\mathbb{R}^{n}$ is continuous and if $\psi$ is
closed, then:

$g_{\psi}(x)=\begin{cases}\max\left\\{\varepsilon\geq
0\,|\,\xi(x,\varepsilon)\vDash\psi\right\\},&\text{ if }\xi(x,0)\vDash\psi,\\\
0,&\text{ if }\xi(x,0)\nvDash\psi.\end{cases}$

###### Proof.

First using Theorem 2.28 of [R+76], the supremum is equal to the maximum if
the set $A=\left\\{\varepsilon\geq
0\,|\,\xi(x,\varepsilon)\vDash\psi\right\\}$ is closed. To show the result, it
is sufficient to show that the set $A$ is closed. Consider a sequence
$\varepsilon_{i}\in\mathbb{R}_{\geq 0}$, $i\in\mathbb{N}$, assume that for all
$i\in\mathbb{N}$, $\xi(x,\varepsilon_{i})\vDash\psi$ and let us show that
$\xi(x,\varepsilon)\vDash\psi$, with
$\varepsilon=\lim_{i\rightarrow\infty}\varepsilon_{i}$. Using the closedness
of the specification $\psi$, one has that
$\lim_{i\rightarrow\infty}\xi(x,\varepsilon_{i})\vDash\psi$, which implies
from the continuity of the map $f$ and using Lemma 7.3 (in Appendix) that
$\lim_{i\rightarrow\infty}\xi(x,\varepsilon_{i})\subseteq\xi(x,\lim_{i\rightarrow\infty}\varepsilon_{i})=\xi(x,\varepsilon)\vDash\psi$.
∎

We also provide sufficient conditions under which the map $g_{\psi}(x)$
remains bounded for an $x\in X$.

######  Proposition 3.4.

Consider the dynamical system $\Sigma$ in (1), an LTLF specification $\psi$
and a point $x\in X$. If the set
$\textsf{L}^{-1}(\psi)\subseteq\mathbb{R}^{n}$ is a bounded subset of
$\mathbb{R}^{n}$, then $g_{\psi}(x)$ is bounded.

###### Proof.

Assume that $g_{\psi}(x)$ is unbounded, hence the set $\xi(x,g_{\psi}(x))$
$\subseteq\textsf{L}^{-1}(\psi)\subseteq\mathbb{R}^{n}$ is unbounded, which
contradicts the fact that the set
$\textsf{L}^{-1}(\psi)\subseteq\mathbb{R}^{n}$ is bounded. Hence,
$g_{\psi}(x)$ is bounded. ∎

###  3.4. Resilience Properties for Linear Systems and Convex Specifications

In this part, we introduce the concept of convex specifications and present an
efficient approach for the computation of the resilience metric for the class
of linear systems and convex specifications.

######  Definition 3.3 (Convex specification).

An LTLF formula $\psi$ is said to be convex if the following holds: for
trajectories $w_{x,i}$, $i\in\\{1,2\\}$, if $w_{x,i}\vDash\psi$, then for any
trajectory $w_{x}=\lambda w_{x,1}+(1-\lambda)w_{x,2}$, $\lambda\in[0,1]$, we
have that $w_{x}\vDash\psi$.

Intuitively, the convexity property states that the specification is preserved
under a convex hull operator. In Appendix 7.3 provides a characterization of
the fragment of LTLF specifications that are convex.

Now, consider the case where the objective is to compute $g_{\psi}(A)$ for a
set $A\subset\mathbb{R}^{n}$. A straightforward approach that was mentioned
earlier in the property (vi) of Proposition 3.1 is to use the fact that
$g_{\psi}(A)=\inf_{x\in A}g_{\psi}(x)$, which requires computing $g_{\psi}(x)$
for all $x\in X$ and that can be computationally infeasible for continuous
sets. In this part, we present an efficient approach to compute $g_{\psi}(A)$
for the case where the set $A$ can be written as a convex closure of a finite
number of points. Our result relies on the superposition principle for linear
systems and the introduced class of convex specifications.

######  Theorem 3.5.

Consider the discrete-time linear system in (1). Consider
$A=\textrm{conv}(c_{1},c_{2},\ldots,c_{p})\subset\mathbb{R}^{n}$ and a convex
specification $\psi$. We have
$g_{\psi}(A)=\min\limits_{i=1,2,\ldots,p}g_{\psi}(c_{i})$.

###### Proof.

First, one has from (vi) in Proposition 3.1 that
$g_{\psi}(X)\leq\min\limits_{i=1,2,\ldots,p}g_{\psi}(c_{i})$. To show that
$g_{\psi}(A)\geq\min\limits_{i=1,2,\ldots,p}g_{\psi}(c_{i})$ , let us denote
$\varepsilon=\min\limits_{i=1,2,\ldots,p}g_{\psi}(c_{i})$ and consider $x\in
A=\textrm{conv}(c_{1},c_{2},\ldots,c_{p})$. We have the existence of
$\alpha_{1},\alpha_{2},\ldots,\alpha_{p}\geq 0$ such that
$\sum\limits_{i=1}^{p}\alpha_{i}c_{i}=x$ and
$\sum\limits_{i=1}^{p}\alpha_{i}=1$. To show the result, let us show that

(7)
$\xi(x,\varepsilon)\subseteq\textrm{conv}\left(\xi(c_{1},\varepsilon),\xi(c_{2},\varepsilon),\ldots,\xi(c_{p},\varepsilon)\right).$

Consider a state trajectory
$w_{x}=x_{0},x_{1},x_{2},\ldots\in\xi(x,\varepsilon)$, we have the existence
of $d_{0},d_{1},d_{2},\ldots\in\Omega_{\varepsilon}(0)$ such that $x_{0}=x$
and $x_{k+1}=Ax_{k}+d_{k}$ for all $k\in\\{0,\ldots,N-1\\}$. Now for
$i\in\\{1,2,\ldots,p\\}$, consider
$w_{x,i}=x_{0,i},x_{1,i},x_{2,i},\ldots\in\xi(c_{i},\varepsilon)$ defined as
follows: $x_{0,i}=c_{i}$ and $x_{k+1,i}=Ax_{k,i}+d_{k}$. Hence, using the fact
that $\sum\limits_{i=1}^{p}\alpha_{i}=1$ one has that for all
$i\in\\{1,2,\ldots,p\\}$ and for all $k\in\\{0,\ldots,N-1\\}$,

$x_{k}=\sum\limits_{i=1}^{p}(\alpha_{i}x_{k,i}+d_{k})\in\textrm{conv}({x_{k,1},x_{k,2},\ldots,x_{k,p}}),$

which in turn implies that
$w_{x}\in\textrm{conv}({w_{x,1},w_{x,2},\ldots,w_{x,p}})$ and (7) holds.

Moreover, since $\varepsilon\leq g_{\psi}(c_{i})$, $i\in\\{1,2,\ldots,p\\}$,
we have from Lemma 3.2 that for all $i\in\\{1,2,\ldots,p\\}$,
$\xi(c_{i},\varepsilon)\vDash\psi$. Hence using the fact that
$\xi(x,\varepsilon)\subseteq\textrm{conv}({\xi(c_{1},\varepsilon),\xi(c_{2},\varepsilon),\ldots,\xi(c_{p},\varepsilon)})$
and from the convexity of the specification $\psi$, we conclude that
$\xi(x,\varepsilon)\vDash\psi$. Hence,
$g_{\psi}(x)=\min\limits_{i=1,2,\ldots,p}g_{\psi}(c_{i})$. ∎

In Sections 4-5, we provide results on the computation of the resilience
metric for linear and nonlinear systems for various specifications.

##  4\. Computation of Resilience for Linear Systems

In this section, we show how the resilience metric can be computed exactly for
some classes of specifications, such as the exact time reachability and
finite-horizon safety. Moreover, we show how it can be approximated, up to any
accuracy, for finite-horizon reachability properties.

###  4.1. Exact-Time Reachability

Consider the linear discrete-time system $\Sigma$ defined by:

(8) $x(k+1)=Ax(k)+d(k),$

with $x(k),d(k)\in\mathbb{R}^{n},k\in\mathbb{N}$ and reachability at a
specific time $N\in\mathbb{N}$ as $\psi=\mathord{\bigcirc}^{N}\Gamma$, for
some polytopic set $\Gamma$. We have the following result showing that the
computation of the resilience metric for linear systems and reachability at a
specific time point specification boils down to a linear optimization problem.

######  Theorem 4.1.

Consider the linear system $\Sigma$ in (8) and the specification
$\psi=\mathord{\bigcirc}^{N}\Gamma$, for $N\in\mathbb{N}$, where $\Gamma$ is a
polytope $\Gamma=\\{x\in X\,|\,Gx\leq H\\}$ with $G\in\mathbb{R}^{q\times n}$
and $H\in\mathbb{R}^{q}$, $q\in\mathbb{N}$. We have

(9) $g_{\psi}(x)\\!=\\!\min\\{\varepsilon\geq 0\,|\,P\geq
0,PA_{b}\\!=\\!E,PB_{b}\leq\varepsilon F(x)\\}$

with

  * •

$A_{b}:=\left[\begin{array}[]{cc}\mathbb{I}\\\
-\mathbb{I}\end{array}\right]\in\mathbb{R}^{2nN\times nN}$ and
$B_{b}:=\left[\begin{array}[]{c}\mathbf{1}\\\
\mathbf{1}\end{array}\right]\in\mathbb{R}^{2nN}$

  * •

$E=[GA^{N-1},GA^{N-2},\ldots,GA,G]\in\mathbb{R}^{q\times nN}$

  * •

$F(x)=H-GA^{N}x\in\mathbb{R}^{q}$.

###### Proof.

We have that

$\displaystyle\begin{cases}g_{\psi}(x)=0\qquad\text{ if }\quad
GA^{N}x>H,\text{ else: }\\\ g_{\psi}(x)=\max\,\\{\varepsilon\geq
0,\,s.t.\,\text{ for all }d_{0},\ldots,d_{N-1}\in\Omega_{\varepsilon}(0)\\\
\qquad\qquad\qquad G(A^{N}x+A^{N-1}d_{0}+\ldots+d_{N-1})\leq H.\\}\end{cases}$

Hence, the computation of the resilience metric $g_{\psi}$ requires solving
the following optimization problem

$\displaystyle\max\\{\varepsilon\geq
0\,|\,EY\leq\frac{1}{\varepsilon}F(x),\text{for all }Y\text{ satisfying
}A_{b}Y\leq B_{b}\\},$

with $Y:=\frac{1}{\varepsilon}[d_{0};d_{1};\ldots,d_{N-1}]$ being the free
variables, $F(x):=H-GA^{N}x$, $E:=[GA^{N-1},GA^{N-2},\ldots,GA,G]$, and
$A_{b}Y\leq B_{b}$ represents the inequality $\|Y\|_{\infty}\leq 1$ in matrix
form using the matrices $A_{b}$ and $B_{b}$ defined above.111In the general
case when the disturbance set has a template in the form of a polytope
$\\{w\,|\,A_{w}w\leq B_{w}\\}$, we have
$A_{b}=\text{diag}(A_{w},A_{w},\ldots,A_{w})$ and
$B_{b}=[B_{w};B_{w};\ldots;B_{w}]$.

$A_{b}=\left[\begin{array}[]{cc}\mathbb{I}\\\
-\mathbb{I}\end{array}\right],B_{b}=\left[\begin{array}[]{c}\mathbf{1}\\\
\mathbf{1}\end{array}\right].$

According to the affine form of Farkas’ Lemma [Sch98], see Theorem 7.2 in
Appendix, one gets the following optimization problem:

$g_{\psi}(x)=\max\\{\varepsilon\geq 0\,|\,P\geq
0,PA_{b}=E,PB_{b}\leq\frac{1}{\varepsilon}F(x)\\},$

which is equivalent to the optimization problem in (9). ∎

We have the following illustrative example for the results of Theorems 4.1 and
3.5.

######  Example 4.1.

Consider the linear dynamical system $\Sigma$ in (8) with
$A=\begin{pmatrix}0.1&-1\\\ -0.5&-0.2\end{pmatrix}$ and the polytope
$\Gamma=\\{x\in X\,|\,Gx\leq H\\}$ with

$G=\begin{pmatrix}0&-0.2747\\\ 0.3714&0\\\ 0&0.3714\\\
-0.2747&0\end{pmatrix},\text{ and }H=\begin{pmatrix}0.9615\\\ 0.9285\\\
0.9285\\\ 0.9615\end{pmatrix}.$

Consider the specification defined by $\psi=\mathord{\bigcirc}^{3}\Gamma$. In
simple words, the objective is to reach the set $\Gamma$ in three time steps,
the set $\Gamma$ is defined as an interval
$\Gamma=[-3.5,3.5]\times[-2.5,2.5]$. Figure 1 represents the resilience metric
$g_{\psi}(x)$ on the domain $A=[-4,6]\times[-4,6]$. Moreover, the numerical
results show that $g_{\psi}(X)=0.1166$. On the other hand, we have that the
set $A=[-4,6]\times[-4,6]=\textrm{conv}(c_{1},c_{2},c_{3},c_{4})$, with
$c_{1}=(-4,-4)$, $c_{2}=(-4,6)$, $c_{3}=(6,-4)$ and $c_{4}=(6,6)$. The
resilience metrics of the vectors $c_{1}$, $c_{2}$, $c_{3}$ and $c_{4}$ are
$g_{\psi}(c_{1})=0.417$, $g_{\psi}(c_{2})=0.1166$, $g_{\psi}(c_{3})=0.1380$
and $g_{\psi}(c_{4})=0.1203$. One can easily check that
$g_{\psi}(X)=\min\limits_{i=1,2,3,4}g_{\psi}(c_{i})=g_{\psi}(c_{2})$, which is
consistent with the result of Theorem 3.5.

Figure 1. Function $g_{\psi}(x)$ for Example 4.1.

###  4.2. Finite-Horizon Safety Specifications

In this section, we provide a closed-form expression of the resilience metric
for the case of linear systems and finite-horizon safety specifications.

######  Theorem 4.2.

Consider the linear system $\Sigma$ in (8) with finite-horizon safety
specification $\psi=\square^{N}\Gamma$, where $\Gamma=\\{x\in X\,|\,Gx\leq
H\\}$, for some $N\in\mathbb{N}$. We have

$\displaystyle g_{\psi}(x)=\max\\{\varepsilon\geq 0\mid P\geq
0,PA_{b}=E,PB_{b}\leq\varepsilon F(x)\\},$

  * •

$A_{b}=\left[\begin{array}[]{cc}\mathbb{I}\\\
-\mathbb{I}\end{array}\right]\in\mathbb{R}^{2nN\times nN}$ and
$B_{b}=\left[\begin{array}[]{c}\mathbf{1}\\\
\mathbf{1}\end{array}\right]\in\mathbb{R}^{2nN}$

  * •

$\displaystyle E=\left[\begin{array}[]{cccccc}G&0&0&\ldots&0&0\\\
GA&G&0&\ldots&0&0\\\ GA^{2}&GA&G&\ldots&0&0\\\
\vdots&\vdots&\vdots&\vdots&\vdots&\vdots\\\
GA^{N-1}&GA^{N-2}&GA^{N-3}&\ldots&GA&G\end{array}\right],$ $\displaystyle
F(x)=[H-GAx\;H-GA^{2}x\;\ldots\;H-GA^{N}x].$

###### Proof.

The resilience metric is given by

$\displaystyle g_{\psi}(x)=$ $\displaystyle\max\,\varepsilon\geq 0,\text{s.t
for all }d_{0},\ldots,d_{N-1}\in\Omega_{\varepsilon}(0),$
$\displaystyle\begin{cases}Gx\leq H\\\ G(Ax+d_{0})\leq H\\\
G(A^{2}x+Ad_{0}+d_{1})\leq H\\\ \vdots\\\
G(A^{N}x+A^{N-1}d_{0}+\ldots+d_{N-1})\leq H.\end{cases}$

Hence, the computation of $g_{\psi}$ requires solving an optimization of the
form

$\displaystyle\max\\{\varepsilon\geq
0\,|\,EY\leq\frac{1}{\varepsilon}F(x),\,\,\forall Y\text{ with }A_{b}Y\leq
B_{b}\\},$

with $Y:=\frac{1}{\varepsilon}[d_{0}^{T},d_{1}^{T},\ldots,d_{N-1}^{T}]^{T}$
being the free variables, where $A_{b}Y\leq B_{b}$ represents the inequality
$\|Y\|_{\infty}\leq 1$ in matrix form. Hence, the result holds using Theorem
7.2 in the Appendix. ∎

Let us remark that by defining $\psi_{i}=\mathord{\bigcirc}^{i}\Gamma$, one
gets $\psi=\wedge_{i=0}^{N}\psi_{i}$. Therefore, it follows from (iii) in
Proposition 3.1 that $g_{\psi}(x)=\min_{i}g_{\psi_{i}}(x)$. Hence, one can
state the previous result in terms of reachability with exact time, with
$g_{\psi_{i}}$ computed previously in Theorem 4.1.

###  4.3. Finite-Horizon Reachability

In this section, we provide two approaches to compute the resilience metric
for linear systems and finite-horizon reachability specifications. The first
approach based on the exact time reachability approach proposed in Section
4.1, and a second approach based on scenario optimization.

####  4.3.1. Computation using LMIs

Consider the linear system $\Sigma$ in (8) with finite-horizon reachability
specification $\psi=\lozenge^{N}\Gamma$ for some set $\Gamma$ as a polytope
$\Gamma:=\\{x\in X\,|\,Gx\leq H\\}$. Then, we have

$\displaystyle g_{\psi}(x)$ $\displaystyle=\max\,\varepsilon\geq
0,\,s.t.\,\text{ for all }d_{0},\ldots,d_{N-1}\in\Omega_{\varepsilon}(0)$ (10)
$\displaystyle\begin{cases}Gx\leq H\text{ or }\\\ G(Ax+d_{0})\leq H\text{ or
}\\\ G(A^{2}x+Ad_{0}+d_{1})\leq H\text{ or }\\\ \vdots\\\
G(A^{N}x+A^{N-1}d_{0}+\ldots+d_{N-1})\leq H.\end{cases}$

in view of (v) in Proposition 3.1, one can select
$\psi_{i}=\mathord{\bigcirc}^{i}\Gamma$, to get $\psi=\vee_{i=0}^{N}\psi_{i}$.
Therefore, we can use the results of reachability in an exact time presented
in Section 4.1 to obtain a lower bound on the resilience metric given by
$g_{\psi}(x)\geq\max_{i}g_{\psi_{i}}(x)$. However, using this approach, we do
not have any measure of how the computed resilience metric is conservative. To
deal with this issue, in the next section, we propose a scenario optimization
approach to compute a lower bound on the resilience metric, together with a
measure of its conservatism.

####  4.3.2. Computation using Scenario Optimization

To get a solution of the optimization problem in (4.3.1). We define
$\mu=\frac{1}{\varepsilon}$, and the set
$\mathcal{Z}=\\{(z_{0},z_{1},\ldots,z_{N-1})\in\Omega_{1}(0)\\}$. Hence, one
gets that $g_{\psi}(x)$ is given by the following optimization problem:

$\displaystyle g_{\psi}(x)=\min\,$ $\displaystyle\mu\geq 0,\,s.t.\,\text{ for
all }z_{0},\ldots,z_{N-1}\in\Omega_{1}(0)$ $\displaystyle\begin{cases}G\mu
x\leq\mu H\text{ or }\\\ G(\mu Ax+z_{0})\leq\mu H\text{ or }\\\ G(\mu
A^{2}x+Az_{0}+z_{1})\leq\mu H\text{ or }\\\ \vdots\\\ G(\mu
A^{N}x+A^{N-1}z_{0}+\ldots+z_{N-1})\leq\mu H.\end{cases}$

If there exists a known $\overline{\mu}\geq 0$ such that the problem is
feasible for $\overline{\mu}$222In this problem, $\bar{\mu}$ can be chosen as
$\bar{\mu}=\frac{1}{\varepsilon}$ for any $\varepsilon>0$ that is solution to
the optimization problem (4.3.1)., the computation of $g_{\psi}(x)$ can be
written in terms of the following non-convex program

$\displaystyle\min\mu$ $\displaystyle\text{ s.t }~{}\mu\in T\text{ and
}\bigvee\limits_{k=1}^{N+1}g_{k}(\mu,z)\leq 0\text{ for all }z\in\mathcal{Z},$

where $T=[0,\overline{\mu}]$ and for $\mu\geq 0$ and
$z=(z_{0},z_{1},\ldots,z_{N-1})\in\mathcal{Z}$, the map $g_{k}$, is defined by
$g_{1}(\mu,z)=G\mu x-\mu H$ and for $k\in\\{2,\ldots,N+1\\}$

(11) $g_{k}(\mu,z)=G(\mu A^{k-1}x+A^{k-2}z_{0}+\ldots+z_{k-2})-\mu H$

The symbol $\bigvee\limits_{k=1}^{N+1}g_{k}(\mu,z)\leq 0$ indicates that there
exists $k\in\\{1,2,\ldots,N\\}$ such that $g_{k}(\mu,z)\leq 0$, for all
$z\in\mathcal{Z}$.

For a chosen parameter $\gamma\in\mathbb{R}^{N+1}_{\geq 0}$, in order to
construct the $\gamma$-scenario program [ESL14] associated with the robust
program defined above, we consider a uniform distribution on the space
$\mathcal{Z}$ and obtain $M$ i.i.d. sample points
$\\{z^{1},z^{2},\ldots,z^{M}\\}$ with
$z^{i}=(z^{i}_{0},z^{i}_{1},\ldots,z^{i}_{N-1})\in\mathcal{Z}$,
$i\in\\{1,2,\ldots,M\\}$. In view of [ESL14], the $\gamma$-scenario program
can be written as

(12) $\displaystyle\min\mu$ $\displaystyle\text{ s.t }~{}\mu\in T\text{ and
}\bigvee\limits_{k=1}^{N}(g_{k}(\mu,z^{i})+\gamma_{k}\leq
0),~{}~{}i\in\\{1,2,\ldots,M\\}.$

We have the following result.

######  Theorem 4.3.

Consider the linear system $\Sigma$ in (8) and the specification
$\psi=\lozenge^{N}\Gamma$, for some $N\in\mathbb{N}$. For $\beta\in(0,1]$ and
$\eta\in[0,1]^{N+1}$, the optimal solution of (12), gives an approximate
optimal solution to the optimization problem (4.3.1) with confidence
$(1-\beta)$, when the number of samples $M$ satisfies

(13) $M\geq
M\left((\eta_{1}^{\frac{1}{N}},\eta_{2}^{\frac{1}{N}},\ldots,\eta_{N+1}^{\frac{1}{N}}),\beta\right)$

and $\gamma=(\gamma_{1},\gamma_{2},\ldots,\gamma_{N+1})$ with
$\gamma_{k}=(\|A\|^{k-1}\|Gx\|+\|H\|)\eta^{\frac{1}{N}}$,
$k\in\\{1,2,\ldots,N+1\\}$, and the map
$M:[0,1]^{N+1}\times(0,1]\rightarrow\mathbb{N}$ is defined for
$\alpha\in[0,1]^{N+1}$ and $\theta\in(0,1]$ by
$M(\alpha,\theta)=\min\\{M\in\mathbb{N}\mid\sum_{k=1}^{N+1}(1-\alpha_{k})^{M}\leq\theta\\}$.
Moreover, the optimal value $J^{*}$ for the optimization problem (4.3.1), and
$J^{*}_{M}$ the optimal value for the optimization problem (12) satisfy
$J^{*}-J^{*}_{N}\leq\|\eta\|_{\infty}$, with confidence $(1-\beta)$.

###### Proof.

We prove the result by using [ESL14, Theorem 4.3]. We first provide a
Lipschitz constant $L_{z,k}$ for the map $(\mu,z)\mapsto g_{k}(\mu,z)$ with
respect to $\mu$, $k\in\\{1,2,\ldots,N+1\\}$. Consider $\mu,\mu^{\prime}\in T$
and $z\in\mathcal{Z}$. We have:

$\displaystyle\|g_{k}(\mu,z)-g_{k}(\mu^{\prime},z)\|$
$\displaystyle\leq\|G(\mu A^{k-1}x+A^{k-2}z_{0}+\ldots+z_{k-2})-\mu H$
$\displaystyle-G(\mu^{\prime}A^{k-1}x+A^{k-2}z_{0}+\ldots+z_{k-2})-\mu^{\prime}H\|$
$\displaystyle+\|(\mu-\mu^{\prime})H\|$
$\displaystyle\leq\|A\|^{k-1}\|(\mu-\mu^{\prime})Gx\|+\|(\mu-\mu^{\prime})H\|$
$\displaystyle\leq(\|A\|^{k-1}\|Gx\|+\|H\|)|\mu-\mu^{\prime}|.$

Hence, $L_{z,k}=\|A\|^{k-1}\|Gx\|+\|H\|$ is a Lipschitz constant for the map
$g_{k}$, $k\in\\{1,2,\ldots,N+1\\}$ with respect to $z$. Moreover, for any
$k\in\\{1,2,\ldots,N+1\\}$, the map$(\mu,z)\mapsto g_{k}(\mu,z)$ is linear and
then convex with respect to $\mu$. Since the distribution on the set
$\mathcal{Z}=\\{(z_{0},z_{1},\ldots,z_{N-1})\in\Omega_{1}(0)\\}$ is uniform,
we choose $h(\varepsilon)=(\frac{\varepsilon}{2})^{N}$. The result follows
from an application [ESL14, Theorem 4.3], [ESL14, Remark 3.9] and [ESL14,
Remark 3.5]. ∎

##  5\. Computation of Resilience for Nonlinear Systems

In this section, we extend the approaches to compute the resilience metric to
the case of nonlinear systems. In particular, we discuss two approaches: (i)
linear optimization-based approach and (ii) SMT-based approach.

###  5.1. Linear Optimization-based Approach

Since the computation of resilience for different types of specifications
relies on reachability with the exact time as a building block, in this
section, we focus on reachability with an exact time specification. The
extension to other specifications can be done following the approaches
presented in the previous section.

Consider the reachability at a specific time point:
$\psi=\mathord{\bigcirc}^{N}\Gamma$, for some set $\Gamma$ defined as a
polytope $\Gamma=\\{x\in X\,|\,Gx\leq H\\}$. We provide a linear optimization-
based solution for computing $g_{\psi}(x)$ for nonlinear systems. we use the
following auxiliary result.

######  Lemma 5.1.

Consider the nonlinear discrete-time system $\Sigma$ defined by:

(14) $x(k+1)=f(x(k)),~{}x(k)\in X\subseteq\mathbb{R}^{n},~{}k\in\mathbb{N}.$

with a continuously differential map $f$. Assume without loss of generality
that $f(0)=0$ and assume the existence of $\overline{\alpha}_{ij}$,
$\underline{\alpha}_{ij}$, $i,j\in\\{1,2,\ldots,n\\}$, such that for all $x\in
X\subset\mathbb{R}^{n}$:

(15) $\underline{\alpha}_{ij}\leq\frac{\partial f_{i}}{\partial
x_{j}}\leq\overline{\alpha}_{ij},~{}i,j\in\\{1,2,\ldots,n\\}$

Then, we have $\underline{A}x\leq f(x)\leq\overline{A}x$ for all $x\in X$,
with $\underline{A},\overline{A}\in\mathbb{R}^{n\times n}$ defined for
$i,j\in\\{1,2,\ldots,n\\}$ by $\underline{A}_{ij}=\underline{\alpha}_{ij}$ and
$\overline{A}_{ij}=\overline{\alpha}_{ij}$.

###### Proof.

Condition (15) implies that for all $x,y\in X$ and for all
$i,j\in\\{1,2,\ldots,n\\}$, we have

$\sum\limits_{j=1}^{n}\underline{\alpha}_{ij}(x_{j}-y_{j})\leq
f_{i}(x)-f_{i}(y)\leq\sum\limits_{j=1}^{n}\overline{\alpha}_{ij}(x_{j}-y_{j}).$

In particular for $y=0$, one gets that for all $x\in X$ and for all
$i\in\\{1,2,\ldots,n\\}$, we have

$\sum\limits_{j=1}^{n}\underline{\alpha}_{ij}x_{j}\leq
f_{i}(x)\leq\sum\limits_{j=1}^{n}\overline{\alpha}_{ij}x_{j}.$

and the result holds. ∎

######  Theorem 5.2.

Consider the nonlinear system $\Sigma$ given by

(16) $x(k+1)=f(x(k))+d(k),~{}x(k)\in
X\subseteq\mathbb{R}^{n},~{}k\in\mathbb{N}$

and the specification $\psi=\mathord{\bigcirc}^{N}\Gamma$, where $\Gamma$ is
the polytope $\Gamma=\\{x\in X\,|\,Gx\leq H\\}$. Assume the existence of
$\overline{\alpha}_{i,j}$, $\underline{\alpha}_{ij}$,
$i,j\in\\{1,2,\ldots,n\\}$, such that for all $x\in X$:

(17) $\underline{\alpha}_{ij}\leq\frac{\partial f_{i}}{\partial
x_{j}}\leq\overline{\alpha}_{ij},~{}i,j\in\\{1,2,\ldots,n\\},$

Consider the matrices $B_{1},B_{2},\ldots,B_{N}\in\mathbb{R}^{n\times n}$
defined as

$B_{k,ij}=\max\left(\sum_{h=1}^{n}G_{ih}\underline{D}_{k,hj},\sum_{h=1}^{c}G_{ih}\overline{D}_{k,hj}\right),$

where $B_{k,{ij}}$ represents the coefficient corresponding to the position
$(i,j)$ of the matrix $B_{k}$, $k\in\\{1,2,\ldots,N\\}$, and the matrices
$\underline{D}_{k},\overline{D}_{k}$ are given with
$\underline{D}_{k}=\underline{A}^{k}$ and $\overline{D}_{k}=\overline{A}^{k}$,
with $\underline{A},\overline{A}$ defined for $i,j\in\\{1,2,\ldots,n\\}$ by
$\underline{A}_{ij}=\underline{\alpha}_{ij}$ and
$\overline{A}_{ij}=\overline{\alpha}_{ij}$. Then, we have

$\displaystyle g_{\psi}(x)\geq\max$ $\displaystyle\\{\varepsilon\geq
0\,|\,B_{N}x+B_{N-1}d_{0}+\ldots+d_{N-1}\leq H$ (18) $\displaystyle\text{ for
all }d(0),\ldots,d(N-1)\in\Omega_{\varepsilon}(0)\\},$

###### Proof.

The computation of resilience metric $g_{\psi}$ for reachability at a specific
time point for the nonlinear system $\Sigma$ requires solving the following
optimization problem

$\displaystyle\max\\{\varepsilon\geq 0\,|\,$ $\displaystyle
Gf(f(\ldots((f(x)+d(0))+d(1)),\ldots)+d(N-1))\leq H,$ (19)
$\displaystyle\text{ for all
}d(0),\ldots,d(N-1)\in\Omega_{\varepsilon}(0)\\}.$

In view of Lemma 5.1, the nonlinear term
$f(f(\ldots(f(x)+d(0))+d(1)),\ldots)+d(N-1))$ can be abstracted into a linear
one as follows:

$\displaystyle\underline{A}^{N}x+\underline{A}^{N-1}d_{0}+\ldots+d_{N-1}\leq$
$\displaystyle f(f(\ldots((f(x)+d(0))+d(1)),\ldots)+d(N-1))\leq$
$\displaystyle\overline{A}^{N}x+\overline{A}^{N-1}d_{0}+\ldots+d_{N-1}.$

It follows that

$\displaystyle Gf(f(\ldots$
$\displaystyle(f(x)+d(0))+d(1)),\ldots)+d(N-1))\leq$ $\displaystyle
B_{N}x+B_{N-1}d_{0}+\ldots+d_{N-1},$

where $B\in\mathbb{R}^{n\times n}$ is defined for $i,j\in\\{1,2,\ldots,n\\}$
by

$B_{k,ij}=\max(\sum\limits_{h=1}^{n}G_{ih}\underline{A}^{k}_{hj},\sum\limits_{h=1}^{c}G_{ih}\overline{A}^{k}_{hj}),$

where $B_{k,{ij}}$ represents the coefficient corresponding to the position
$(i,j)$ of the matrix $B_{k}$, $k\in\\{1,2,\ldots,N\\}$. Hence, the
optimization problem in (5.1) can be relaxed into the following linear
optimization problem.

$\displaystyle\max\\{\varepsilon\geq 0\,|\,$ $\displaystyle
B_{N}x+B_{N-1}d_{0}+\ldots+d_{N-1}\leq H$ (20) $\displaystyle\text{ for all
}d(0),\ldots,d(N-1)\in\Omega_{\varepsilon}(0)\\},$

which ends the proof. ∎

Theorem 5.2 shows how to transform the problem of computing the resilience
metric for nonlinear systems and exact-time reachability, into a problem for a
linear system that can be resolved using the approach proposed in Section 4.1.

###  5.2. SMT-based Approach

This section provides an SMT-based solution to compute resilience for
nonlinear systems without any linear approximation. In addition, the approach
also relaxes polytopic assumptions over subsets of state space corresponding
to atomic propositions. To utilize an SMT-based solution, we consider the
complement of the resilience definition (3) as given below. For a given LTLF
formula $\psi$, the resilience in (3) can be equivalently written as

(21) $g_{\psi}(x)=\begin{cases}\inf\left\\{\varepsilon\geq
0\,|\,\xi(x,\varepsilon)\nvDash\psi\right\\},&\text{ if
}\xi(x,0)\nvDash\psi\\\ 0,&\text{ if }\xi(x,0)\vDash\psi,\end{cases}$

where $\neg\psi$ is a negation of the formula $\psi$.
Consider the set $\Gamma$ represented by logical conjunctions and disjunctions
of $P$ nonlinear predicate functions $p_{i}:X\rightarrow\mathbb{R}$ as given
below:

$\displaystyle\Gamma=\\{x\in X\mid p_{i}(x)\leq 0\ \Theta_{i}\ p_{i+1}(x)\leq
0\text{ is {true}},$ (22) $\displaystyle\forall
i\in\\{1,\ldots,P\\},\Theta_{i}\in\\{\vee,\wedge\\}\\}.$

Now, consider LTLF formula $\psi=\mathord{\bigcirc}^{N}\Gamma$, one can
compute the resilience $g_{\psi}(x)$ in (21) by solving the following
optimization problem:

$g_{\psi}(x)=\min\,\varepsilon\geq 0,\,s.t.\,\text{ there exists
}d_{0},\ldots,d_{N-1}\in\Omega_{\varepsilon}(0)\text{ satisfying
}\neg\big{(}p_{i}(\Psi)\leq 0\Theta_{i}p_{i+1}(\Psi)\leq 0\big{)}\text{ is
{true}},\forall
i\in\\{1,\ldots,P\\},\quad\quad\Theta_{i}\in\\{\vee,\wedge\\}$, where
$\Psi=f(f(\ldots((f(x)+d(0))+d(1)),\ldots)+d(N-1))$.

For a given value of $\varepsilon$, one can leverage SMT solvers such as dReal
[GKC13] to check the feasibility of the condition. Then, the near minimum
value of $\varepsilon$ can be found using the bisection method. In a similar
way, one can formulate optimization problems for other formulas
$\psi=\lozenge^{N}\Gamma$ and $\psi=\square^{N}\Gamma$.

##  6\. Case Studies

Temperature Regulation. We consider the problem of regulating the temperature
in a circular building of $9$ rooms. The dynamics of the room temperatures are
given by

$\displaystyle{T}_{i}(k+1)$
$\displaystyle=T_{i}(k)+\alpha(T_{i+1}(k)+T_{i-1}(k)-2T_{i}(k))$
$\displaystyle+\beta(T_{e}+\delta T_{e}-T_{i}(k)),\qquad
i\in\\{1,2,\ldots,9\\},$

where $T_{i+1}$ and $T_{i-1}$ are the temperatures of the neighbor rooms (here
$T_{0}=T_{9}$ and $T_{9+1}=T_{1}$), $T_{e}=0^{\circ}C$ is the outside
temperature, considered as a disturbance and $\alpha$ and $\beta$ are the
conduction factors. The numerical parameters are taken from [SGF20] and given
by $\alpha=0.45$ and $\beta=0.045$.

We consider two scenarios for this example. In the first scenario, the desired
behavior of the system is as follows: The temperatures of the $9$ rooms
initiated in the set $X_{0}=[24,25]^{9}$ should reach the target set
$X_{T}=[21,22]^{9}$ exactly at $N=3$ steps while remaining in the safe set
$X_{S}=[20.5,25]^{9}$. This behavior can be described by the LTLF formula:

(23) $\psi=\psi_{1}\wedge\psi_{2},\text{ with }\psi_{1}=\square^{3}X_{S}\text{
and }\psi_{2}=\mathord{\bigcirc}^{3}X_{T}.$

The objective is to compute the range of admissible external disturbances
$\delta T_{e}$ under which any trajectory of the system initiated in the set
$X_{0}$ satisfies $\psi$. Since the system is linear and the set of initial
states $X_{0}$ is convex, we rely on Theorems 4.1, 4.2 and 3.5 to compute the
resilience metric. The numerical implementations show that the resilience
metric is given by the set $\Delta_{T_{e}}=[-2.23,2.23]^{\circ}C$.

In the second scenario, the desired behavior of the system is as follows: the
temperatures of the $9$ rooms initiated in the set $X_{0}=[50,51]^{9}$ should
remain in the safe set $X_{S}=[10,51]^{9}$, reach the target set
$X_{T_{1}}=[41,43]^{9}$ exactly at $N=4$ steps. Moreover, once the set
$X_{T_{2}}=[31,35]^{9}$ is reached, the temprature should reach the target set
$X_{T_{3}}=[14,28]^{9}$ exactly at $N=11$ steps. This behavior can be
described by the LTLF formula: $\psi=\psi_{1}\wedge\psi_{2}\wedge\psi_{3}$,
with

(24) $\psi_{1}=\square^{2}X_{S},\psi_{2}=\mathord{\bigcirc}^{4}X_{T_{1}}\text{
and }\psi_{3}=X_{T_{2}}\implies\mathord{\bigcirc}^{11}X_{T_{3}}.$

The objective is to compute the range of admissible external disturbances
$\delta T_{e}$ under which any trajectory of the system initiated in the set
$X_{0}$ satisfies $\psi$. The numerical implementations show that the
resilience metric is given by the set $\Delta_{T_{e}}=[-3.12,3.12]^{\circ}C$.

Figure 2 (top) shows the nominal trajectories (with $\delta T_{e}=0$) of the
system for the $9$ rooms. Figure 2 (bottom) shows the trajectories of the
system for the $9$ rooms with a disturbance $\delta T_{e}$ randomly chosen in
the set $\Delta_{T_{e}}=[-3.12,3.12]^{\circ}C$.

Figure 2. Evolution of the temperatures in the nine rooms with a disturbance
$\delta T_{e}=0$ (top) and with a disturbance $\delta T_{e}$ randomly chosen
in the set $\Delta_{T_{e}}=[-3.12,3.12]^{\circ}C$ (bottom). The red boundaries
represent the safe set $X_{S}$. The green region represents the set
$X_{T_{1}}$. The black region represents the set $X_{T_{2}}$ and the blue
region represents the set $X_{T_{3}}$.

Adaptive Cruise Control. Consider a vehicle moving along a straight road. The
dynamics of the vehicle is adapted from [SGF20] and given by the following
difference equation:

(25) $v(k+1)=v(k)+\frac{\tau}{m}(f_{0}-f_{1}v(k)-f_{2}v(k)^{2}),$

where $m>0$ is the mass of the vehicle, $v\geq 0$ represents the velocity of
the vehicle and the term $f_{0}-f_{1}-f_{2}v^{2}$ includes the rolling
resistance and aerodynamics and $\tau$ represents a sampling period. Moreover,
we include a lead vehicle whose velocity $v_{0}\geq 0$ is constant. The
dynamics of the system are

$\left\\{\begin{array}[]{r c l}h(k+1)&=&h(k)+\tau(v_{0}+\delta v_{0}-v(k))\\\
v(k+1)&=&v(k)+\frac{\tau}{m}(f_{0}+\delta
f_{0}-f_{1}v(k)-f_{2}v(k)^{2}),\end{array}\right.$

where $\delta v_{0}$ is the uncertainty on the velocity $v_{0}$ of the lead
vehicle and $\delta f_{0}$ is the uncertainty on the parameter $f_{0}$. The
desired behavior can be described by the following LTLF formula:
$\psi:=\psi_{1}\wedge\psi_{2}$ with $\psi_{1}:=\square^{4}X_{S}$ and
$\psi_{2}:=\mathord{\bigcirc}^{4}X_{T_{1}}$. This behavior can be interpreted
as follows: the relative position should remain in the safe set
$X_{S}=[h_{S,\min},h_{S,\max}]\times[v_{S,\min},v_{S,\max}]$ and the relative
position and velocity should reach the set
$X_{T_{1}}=[h_{T_{1},\min},h_{T_{1},\max}]\times[v_{T_{1},\min},v_{T_{1},\max}]$
in $4$ steps.

The objective is to compute the resilience metric under which the trajectory
of the system initiated from $x_{0}=(60,15)$ satisfies $\psi$. The numerical
values of the vehicle parameters and the considered specifications are given
in Table 1.

We use the approaches developed in Section 5 to deal with this nonlinear
system. First, we used linear optimization based approach. One can easily
check that the values of the parameters $\alpha_{i,j}$, $i,j\in\\{1,2\\}$,
given by $\underline{\alpha}_{11}=\overline{\alpha}_{11}=1,$
$\underline{\alpha}_{12}=\overline{\alpha}_{12}=-1$,
$\underline{\alpha}_{21}=\overline{\alpha}_{21}=0$,
$\underline{\alpha}_{22}=1-\frac{\tau}{m}-2f_{2}v_{S,\max}$ and
$\overline{\alpha}_{22}=1-\frac{\tau}{m}-2f_{2}v_{S,\min}$ satisfy the
inequalities in (17), where the bounds on $\underline{\alpha}_{22}$ and
$\overline{\alpha}_{22}$ follows from the fact that we are interested in
dealing with velocities of the following vehicle within the interval
$[v_{S,\min},v_{S,\max}]$. Hence, in view of Theorem 5.2, one can use the
linear program in (5.2) to compute an approximation of the resilience metric
$g_{\psi}$. The numerical implementations show that the value of the
resilience metric is given by $g_{\psi}=0.0042$. Next, we computed resilience
$g_{\psi}$ using SMT solver dReal and a bisection approach as discussed in
Section 5.2. The obtained value of $g_{\psi}$ using SMT based approach is
$0.0121$ which is less conservative compared to that of using linearization
based approach. The corresponding values for the admissible resilience metric
on the variable $\delta v_{0}$ and $\delta f_{0}$ are given by, $\Delta
v_{0}=[-0.012,0.012]$ and $\Delta f_{0}=[-16.6,16.6]$.

###  6.1. DC Motor

Our third case study is a $2$-dimensional discrete-time DC motor adapted from
[Ade13] as follows:

$\displaystyle x_{1}(k+1)$
$\displaystyle=x_{1}(k)+\tau\big{(}\frac{-R}{L}x_{1}(k)-\frac{k_{dc}}{L}x_{2}(k)\big{)},$
$\displaystyle x_{2}(k+1)$
$\displaystyle=x_{2}(k)+\tau\big{(}\frac{-k_{dc}}{J}x_{1}(k)-\frac{b}{J}x_{2}(k)+d(k)\big{)},$

where $x_{1},x_{2},R=1,L=0.5$, and $J=0.01$ are the armature current, the
rotational speed of the shaft, the electric resistance, the electric
inductance, and the moment of inertia of the rotor, respectively. In addition,
$\tau=0.01,b=0.1,$ and $K_{dc}=0.01$ represent both the motor torque and the
back electromotive force. The desired behavior of the system is as follows:
the armature current $(x_{1})$ and the rotational speed of the shaft $(x_{2})$
starting from initial state $x_{0}=(0.4,0.5)$ should reach the target set
$X_{T}=[0.1,0.35]\times[0.1,0.35]$ within $3$ steps. The corresponding LTLF
formula is $\psi=\lozenge^{3}X_{T}$. The objective is to compute the maximal
value of the disturbance $d$ under which any trajectory of the system
satisfies the specification $\psi$.

We use the Scenario approach proposed in Section 4.3.2, with
$\eta=(0.01,0.01,0.01,0.01)$ and $\beta=0.01$. To achieve the chosen
confidence and accuracy, we choose a sampling size $M$ according to Theorem
4.3. The numerical implementations show that, with confidence $(1-\beta)$, the
resilience metric is given by $g_{\psi}(x_{0})=0.0485$.

Figure 3 shows examples of $1000$ trajectories under disturbance trajectories
chosen randomly on the interval $[-g_{\psi}(x_{0}),g_{\psi}(x_{0})]^{2}$. In
green, we represent the boundaries of the target set $X_{T}$. One can readily
see that almost all the trajectories satisfy the specification $\psi$. Since
the obtained value of the resilience metric $g_{\psi}(x_{0})$ is valid with
confidence $(1-\beta)$, one can also see that some trajectories fail to
satisfy the specification $\psi$.

Figure 3. Evolution of $1000$ perturbed trajectory starting from the initial
condition $x_{0}=(0.4,0.5)$. The green boundaries represent the target set
$X_{T}$.

##  7\. Conclusions and Discussion

We provided a new resilience metric for cyber-physical systems that integrates
the dynamical evolution of the system with temporal logic requirements. We
showed how this resilience metric can be computed for linear and nonlinear
discrete-time models of the system and instances of linear temporal logic
specifications. The related computations utilize Farkas’ Lemma, linear
programs, scenario optimization and SMT-based approaches. A more detailed
analysis is provided for the new introduced classes of closed and convex
finite linear temporal logic specifications. Illustrating our findings, we
showcase computational demonstrations spanning temperature regulation in
buildings, adaptive cruise control, and DC motors.

In the future, we plan to develop techniques for enhancing resilience (i.e.,
maximizing resilience over decision variables) and extending the ideas to
continuous-time dynamical systems. We also plan to use the proposed resilience
metric for the design of contracts for large-scale networked systems [SGF21]
by specifying the largest class of disturbances that a subsystem in the
network can tolerate.

Our new notion of resilience can be used by other research communities working
on the safety and security of critical infrastructures and safety-critical
systems. For instance, it will be useful for designing minimally-invasive
load-altering attacks [KP21] and providing a vulnerability map with respect to
logical requirements on the behavior of the power network [WVS20].

Table 1. Vehicle and safety parameters Parameter | Value | Unit
---|---|---
$M$ | $1370$ | $Kg$
$f_{0}$ | $51.0709$ | $N$
$f_{1}$ | $5$ | $Ns/m$
$f_{2}$ | $0.4161$ | $Ns^{2}/{m^{2}}$
$h_{T_{1},\min}$ | $60.2$ | $m$
$h_{T_{1},\max}$ | $60.5$ | $m$
$v_{T_{1},\min}$ | $14.6$ | $m/s$
$v_{T_{1},\max}$ | $14.9$ | $m/s$
$h_{S,\min}$ | $59.7$ | $m$
$h_{S,\max}$ | $60.7$ | $m$
$v_{S,\min}$ | $14.5$ | $m/s$
$v_{S,\max}$ | $15.5$ | $m/s$

## APPENDIX

###  7.1. Auxilliary results

The following theorem is borrowed from the book [Sch98, Corollary 7.lh] and is
known as the affine form of Farkas’ lemma.

######  Theorem 7.1 ([Sch98]).

Let the system $Ax\leq b$ of linear inequalities have at least one solution,
and suppose that the linear inequality $cx\leq\delta$ holds for each $x$
satisfying $Ax\leq b$. Then for some $\delta^{\prime}\leq\delta$ the linear
inequality $cx\leq\delta^{\prime}$ is a nonnegative linear combination of the
inequalities in the system $Ax\leq b$.

######  Theorem 7.2.

Suppose the set $Ex\leq F$ is not empty. The following two statements are
equivalent:

  * •

$Ex\leq F$ holds for all $x$ with $Ax\leq b$;

  * •

There exists a non-negative matrix $P$ such that $PA=E$ and $Pb\leq F$.

###### Proof.

We use the above theorem for each row of $Ex\leq F$ indicated by $E_{i}x\leq
F_{i}$. Then there is nonnegative row vectors $\lambda_{i}\geq 0$ such that
$\lambda_{i}A=E_{i}$ and $\lambda_{i}b\leq F_{i}$. The non-negative constants
in $\lambda_{i}$ are referred to as the Farkas multipliers. Putting all these
multipliers as row vectors in a matrix $P$ we get that $PA=E$ with $Pb\leq F$
∎

######  Lemma 7.3.

Consider the discrete-time system $\Sigma$ in (1). If the map
$f:\mathbb{R}^{n}\rightarrow\mathbb{R}^{n}$ is continuous, then for any
$x\in\mathbb{R}^{n}$, and for any sequence $\varepsilon_{i}\geq 0$,
$i\in\mathbb{N}$, the map $\xi(x,.):\mathbb{R}_{\geq
0}\rightrightarrows\mathbb{R}^{n}$, defined in (2) satisfies:

(26)
$\lim_{i\rightarrow\infty}\xi(x,\varepsilon_{i})\subseteq\xi(x,\lim_{i\rightarrow\infty}\varepsilon_{i})$

###### Proof.

The proof follows easily from the continuity of the map $f$ and the fact that
$\lim_{i\rightarrow\infty}\Omega_{\varepsilon_{i}}(0)=\Omega_{\lim_{i\rightarrow\infty}\varepsilon_{i}}(0)$
∎

The condition in (26) can be seen as a sufficient condition for the outer
semicontinuity [Aub09] of the map $\xi$ with respect to its second variable.

###  7.2. Closed specifications

In this part, we provide some characterizations of the concept of closed
specifications introduced in Definition 3.2. We first relate the closedness of
the atomic proposition $p$ to the closedness of the corresponding state space
$L^{-1}(p)\subseteq X$.

######  Proposition 7.4.

Consider the alphabet $\Sigma_{a}=2^{AP}$ with the labeling function
$L:X\rightarrow\Sigma_{a}$. Then the atomic proposition $p\in AP$ is closed if
the set $L^{-1}(p):=\cup\\{L^{-1}(A)|A\in\Sigma_{a}\text{ such that }p\in
A\\}$ is a closed subset of $X$. Moreover, the specification $\neg p$, with
$p$ being the atomic proposition, is closed if the set $L^{-1}(p)$ is an open
subset of $X$.

###### Proof.

Consider a converging sequence of trajectories
$w_{x,i}=(x_{0,i},x_{1,i},\ldots,x_{N-1,i})$, $i\in\mathbb{N}$, such that for
all $i\in\mathbb{N}$, $w_{x,i}\vDash\psi$. Hence, one gets that for all
$i\in\mathbb{N}$, $x_{0,i}\in L^{-1}(\psi)$. It follows from the closedness of
the set $L^{-1}(\psi)$, that $x_{0}=\lim_{i\rightarrow\infty}x_{0,i}\in
L^{-1}(\psi)$ which in turn implies that
$w_{x}=\lim_{i\rightarrow\infty}w_{x,i}\vDash\psi$ and that $\psi$ is a closed
specification. The proof of the second result follows immediately. ∎

Now, we provide a characterization of the fragement of LTLf specifications
that are closed.

######  Proposition 7.5.

Let $\psi_{1},\psi_{2}$ be closed specifications and $M\in\mathbb{N}$, with
$M\leq N$, where $N$ represents the length of the trajectories of the LTLF
specifications $\psi_{1},\psi_{2}$. Let $\psi$ be constructed with the grammar

$\displaystyle\psi:=$
$\displaystyle\psi_{1}\,|\,\psi_{1}\wedge\psi_{2}\,|\,\psi_{1}\vee\psi_{2}\,|\,\mathord{\bigcirc}\psi_{1}\,|\,\mathord{\bigcirc}^{M}\psi_{1}\,|\,\square\psi_{1}\,|\,\square^{M}\psi_{1}\,|\,\lozenge\psi_{1}\,|\,\lozenge^{M}\psi_{1}.$

Then, the specification $\psi$ is a closed specification.

###### Proof.

We provide proof for each property separately.

Proof for $\psi_{1}\wedge\psi_{2}$: Consider a converging sequence of
trajectories $w_{x,i}$, $i\in\mathbb{N}$, such that for all $i\in\mathbb{N}$,
$w_{x,i}\vDash\psi_{1}$ and $w_{x,i}\vDash\psi_{2}$, then from closedness of
the specifications $\psi_{1}$ and $\psi_{2}$, it follows that
$w_{x}=\lim_{i\rightarrow\infty}w_{x,i}\vDash\psi_{1}$ and
$w_{x}\vDash\psi_{2}$, which implies that $w_{x}\vDash\psi_{1}\wedge\psi_{2}$.

Proof for $\psi_{1}\vee\psi_{2}$: Consider a converging sequence of
trajectories $w_{x,i}$, $i\in\mathbb{N}$, such that for all $i\in\mathbb{N}$,
$w_{x,i}\vDash\psi_{1}$ or $w_{x,i}\vDash\psi_{2}$. Thus, we can use this fact
to construct two subsequences: the subsequence $w_{x,j_{k}}$ where
$w_{x,j_{k}}$ is a term of $w_{x,i}$ satisfying $w_{x,j_{k}}\vDash\psi_{1}$,
and the subsequence $w_{x,h_{k}}$ where $w_{x,h_{k}}$ is a term of $w_{x,i}$
satisfying $w_{x,h_{k}}\vDash\psi_{2}$. First, since all the terms of the
subsequence $w_{x,j_{k}}$ satisfy $w_{x,j_{k}}\vDash\psi_{1}$, one gets from
the closedeness of $\psi_{1}$ that
$w^{1}=\lim_{k\rightarrow\infty}w_{x,j_{k}}\vDash\psi_{1}$. Similarly, one
gets that $w^{2}=\lim_{k\rightarrow\infty}w_{x,h_{k}}\vDash\psi_{2}$. Now,
since the sequence $w_{x,i}$, $i\in\mathbb{N}$, is converging and has one
limit, the limit of $w_{x,i}$ is either $w^{1}$ or $w^{2}$. Hence, one gets
that $w_{x}=\lim_{i\rightarrow\infty}\vDash\psi_{1}\vee\psi_{2}$.

Proof for $\mathord{\bigcirc}^{M}\psi_{1}$: Consider a converging sequence of
trajectories $w_{x,i}$, $i\in\mathbb{N}$, such that for all $i\in\mathbb{N}$,
$w_{x,i}\vDash\mathord{\bigcirc}^{M}\psi_{1}$, with
$w_{x,i}=(x_{0,i},x_{1,i},\ldots,x_{N-1,i})$. Hence, one gets that for all
$i\in\mathbb{N}$, $x_{M,i}\vDash\psi_{1}$, using the fact that
$w_{x}=\lim_{i\rightarrow\infty}w_{x,i}=x_{0},x_{1},x_{2},\ldots,x_{N}$, and
that $x_{M}=\lim_{i\rightarrow\infty}x_{M,i}$, once gets from the closedness
of the specification $\psi_{1}$ that
$x_{M}=\lim_{i\rightarrow\infty}x_{M,i}\vDash\psi_{1}$, which in turn implies
that $w_{x}\vDash\mathord{\bigcirc}^{M}\psi_{1}$. The proof for the case of
$\mathord{\bigcirc}\psi_{1}$ can be derived similarly.

Proof for $\square^{M}\psi_{1}$: Consider a converging sequence of
trajectories $w_{x,i}$, $i\in\mathbb{N}$, such that for all $i\in\mathbb{N}$,
$w_{x,i}\vDash\square^{M}\psi_{1}$, with
$w_{x,i}=(x_{0,i},x_{1,i},\ldots,x_{N-1,i})$. Hence, one gets that for all
$i\in\mathbb{N}$ and for all $j\in\\{0,1,\ldots M-1\\}$,
$x_{j,i}\vDash\psi_{1}$. Using the fact that
$w_{x}=\lim_{i\rightarrow\infty}w_{x,i}=x_{0},x_{1},x_{2},\ldots,x_{N}$, and
that $x_{j}=\lim_{i\rightarrow\infty}x_{j,i}$, $j\in\\{0,1,\ldots,M-1\\}$,
once gets from the closedness of the specification $\psi_{1}$ that
$x_{j}=\lim_{i\rightarrow\infty}x_{j,i}\vDash\psi_{1}$,
$j\in\\{0,1,\ldots,M-1\\}$, which in turn implies that
$w_{x}\vDash\square^{M}\psi_{1}$. The proof for the case of $\square\psi_{1}$
can be derived similarly.

Proof for $\lozenge^{M}\psi_{1}$: Consider a converging sequence of
trajectories $w_{x,i}$, $i\in\mathbb{N}$, such that for all $i\in\mathbb{N}$,
$w_{x,i}\vDash\lozenge^{M}\psi_{1}$. Hence, for all $i\in\mathbb{N}$, we have
the existence of $j\in\\{0,1,\ldots,M-1\\}$ such that
$w_{x,i}\vDash\mathord{\bigcirc}^{j}\psi_{1}$. Thus, we can use this fact to
construct $M$ subsequence, where the kth subsequence $w_{x,k_{h}}$ is made of
the terms of $w_{x,i}$ satisfying
$w_{x,k_{h}}\vDash\mathord{\bigcirc}^{k}\psi_{1}$, $k\in\\{0,1,\ldots,M-1\\}$.
First, since all the terms of the subsequence $w_{x,k_{h}}$ satisfy
$w_{x,k_{h}}\vDash\mathord{\bigcirc}^{k}\psi_{1}$, one gets from the
closedenss of the specification $\psi_{1}$ and the fact that the closedeness
property is preserved under $\mathord{\bigcirc}^{k}$ that
$w^{k}=\lim_{h\rightarrow\infty}w_{x,k_{h}}\vDash\mathord{\bigcirc}^{k}\psi_{1}$.
Now, since the sequence $w_{x,i}$, $i\in\mathbb{N}$, is converging and has one
limit, we have the existence of $k\in\\{0,1,\ldots,M-1\\}$ such that
$w_{x}=\lim_{i\rightarrow\infty}=w^{k}$. Using the fact that
$w^{k}\vDash\mathord{\bigcirc}^{k}\psi_{1}$ for all $k\in\\{0,1,\ldots,M-1\\}$
one gets that $w_{x}\vDash\lozenge^{M}\psi_{1}$. ∎

###  7.3. Convex specifications

In this part, we provide some characterizations of the concept of convex
specifications introduced in Definition 3.3. We first relate the convexity of
the atomic proposition $p$ to the convexity of the corresponding state space
$L^{-1}(p)\subseteq X$.

######  Proposition 7.6.

Consider the alphabet $\Sigma_{a}=2^{AP}$ with the labeling function
$L:X\rightarrow\Sigma_{a}$. Then the atomic proposition $p\in AP$ is convex if
the set $L^{-1}(p):=\cup\\{L^{-1}(A)|A\in\Sigma_{a}\text{ such that }p\in
A\\}$ is a convex set.

###### Proof.

Consider the trajectories $w_{x,i}$, $i\in\\{1,2\\}$, such that
$w_{x,i}\vDash\mathord{\bigcirc}^{M}\psi$ and consider $\lambda\in[0,1]$.
Define $w_{x}=\lambda w_{x,1}+(1-\lambda)w_{x,2}$. For
$w_{x,i}=(x_{0,i},x_{1,i},\ldots,x_{N-1,i})$, $i=1,2$, one gets that
$x_{0,i}\vDash\psi$. It follows from the convexity of the set $L^{-1}(\psi)$,
that $x_{0}=\lambda x_{0,1}+(1-\lambda)x_{0,2}$ which in turn implies that
$w_{x}=\lambda w_{x,1}+(1-\lambda)w_{x,2}\vDash\psi$ and that $\psi$ is a
convex specification. ∎

Now, we provide a characterization of the fragment of LTLf specifications that
are convex.

######  Proposition 7.7.

Let $\psi_{1},\psi_{2}$ be convex specification and $M\in\mathbb{N}$, with
$M\leq N$, where $N$ represents the length of the trajectories of the LTLF
specifications $\psi_{1},\psi_{2}$. Let $\psi$ be constructed with the grammar

$\displaystyle\psi:=\psi_{1}\,|\,\psi_{1}\wedge\psi_{2}\,|\,\mathord{\bigcirc}\psi_{1}\,|\,\mathord{\bigcirc}^{M}\psi_{1}\,|\,\square\psi_{1}\,|\,\square^{M}\psi_{1}.$

Then, the specification $\psi$ is a convex specification.

###### Proof.

Proof for $\psi_{1}\wedge\psi_{2}$: Consider the trajectories $w_{x,i}$,
$i\in\\{1,2\\}$, such that $w_{x,i}\vDash\psi_{1}\wedge\psi_{2}$ and consider
$\lambda\in[0,1]$. Define $w_{x}=\lambda w_{x,1}+(1-\lambda)w_{x,2}$. Since
both $\psi_{1}$ and $\psi_{2}$ are convex, we have that $w_{x}\vDash\psi_{1}$
and $w_{x}\vDash\psi_{2}$ and the result holds.

Proof for $\mathord{\bigcirc}^{M}\psi_{1}$: Consider the trajectories
$w_{x,i}$, $i\in\\{1,2\\}$, such that
$w_{x,i}\vDash\mathord{\bigcirc}^{M}\psi_{1}$ and consider $\lambda\in[0,1]$.
Define $w_{x}=\lambda w_{x,1}+(1-\lambda)w_{x,2}$. For
$w_{x,i}=(x_{0,i},x_{1,i},\ldots,x_{N-1,i})$, $i=1,2$, one gets that
$x_{M,i}\vDash\psi_{1}$. Since $\psi_{1}$ is a convex specification, it
follows that $w_{x,N}=\lambda x_{M,1}+(1-\lambda)x_{M,2}$, with
$w_{x}=x_{1},x_{2},\ldots,x_{N}$, which in turn implies that
$w_{x}\vDash\mathord{\bigcirc}^{M}\psi_{1}$. The proof for the case of
$\mathord{\bigcirc}\psi_{1}$ can be derived similarly.

Proof for $\square^{M}\psi_{1}$: Consider the trajectories $w_{x,i}$,
$i\in\\{1,2\\}$, with $w_{x,i}=(x_{0,i},x_{1,i},\ldots,x_{N-1,i})$ and such
that $w_{x,i}\vDash\square^{M}\psi_{1}$ and consider $\lambda\in[0,1]$. Define
$w_{x}=\lambda w_{x,1}+(1-\lambda)w_{x,2}$. Hence, one gets that for $i=1,2$
and for all $j\in\\{0,1,\ldots M-1\\}$, $x_{j,i}\vDash\psi_{1}$. Using the
fact that $\psi_{1}$ is a convex specification, it follows that the trajectory
$w_{x}=(x_{0},x_{1},\ldots,x_{N-1})$ satisfies for $j\in\\{0,1,\ldots M-1\\}$
that $x_{j}=\lambda x_{j,1}+(1-\lambda)x_{j,2}\vDash\psi_{1}$, which in turn
implies that $w_{x}\vDash\square^{M}\psi_{1}$. The proof for the case of
$\square\psi_{1}$ can be derived similarly. ∎

## References

  * [ABK14] Shaull Almagor, Udi Boker, and Orna Kupferman.  Discounting in LTL.  In International Conference on Tools and Algorithms for the Construction and Analysis of Systems, pages 424–439. Springer, 2014.
  * [Ade13] Philip. A. Adewuyi.  DC motor speed control: A case between PID controller and fuzzy logic controller.  international journal of multidisciplinary sciences and engineering, 4(4):36–40, 2013.
  * [ATMMR18] Arend Aerts, Bryan Tong Minh, Mohammad Reza Mousavi, and Michel A. Reniers.  Temporal logic falsification of cyber-physical systems: An input-signal-space optimization approach.  In 2018 IEEE International Conference on Software Testing, Verification and Validation Workshops (ICSTW), pages 214–223. IEEE, 2018.
  * [Aub09] J.-P. Aubin.  Viability theory.  Springer Science & Business Media, 2009.
  * [BDD+18] Ezio Bartocci, Jyotirmoy Deshmukh, Alexandre Donzé, Georgios Fainekos, Oded Maler, Dejan Ničković, and Sriram Sankaranarayanan.  Specification-based monitoring of cyber-physical systems: a survey on theory, tools and applications.  In Lectures on Runtime Verification, pages 135–175. 2018.
  * [BK08] Christel Baier and Joost-Pieter Katoen.  Principles of model checking.  MIT press, 2008.
  * [BMM14] Patricia Bouyer, Nicolas Markey, and Raj Mohan Matteplackel.  Averaging in LTL.  In Concurrency Theory, pages 266–280. Springer, 2014.
  * [BT18] Clark Barrett and Cesare Tinelli.  Satisfiability modulo theories.  Springer, 2018.
  * [CLSP22] Hongkai Chen, Shan Lin, Scott A Smolka, and Nicola Paoletti.  An STL-based formulation of resilience in cyber-physical systems.  In International Conference on Formal Modeling and Analysis of Timed Systems, pages 117–135. Springer, 2022.
  * [CSPL23] Hongkai Chen, Scott A Smolka, Nicola Paoletti, and Shan Lin.  An STL-based approach to resilient control for cyber-physical systems.  In Hybrid Systems: Computation and Control (HSCC), pages 1–12, 2023.
  * [DFM13] Alexandre Donzé, Thomas Ferrere, and Oded Maler.  Efficient robust monitoring for stl.  In Computer Aided Verification: 25th International Conference (CAV), pages 264–279. Springer, 2013.
  * [DM10] Alexandre Donzé and Oded Maler.  Robust satisfaction of temporal logic over real-valued signals.  In International Conference on Formal Modeling and Analysis of Timed Systems, pages 92–106. Springer, 2010.
  * [ESL14] Peyman Mohajerin Esfahani, Tobias Sutter, and John Lygeros.  Performance bounds for the scenario approach and an extension to a class of non-convex programs.  IEEE Transactions on Automatic Control, 60(1):46–58, 2014.
  * [FMPS18] Samira S Farahani, Rupak Majumdar, Vinayak S Prabhu, and Sadegh Soudjani.  Shrinking horizon model predictive control with signal temporal logic constraints under stochastic disturbances.  IEEE Transactions on Automatic Control, 64(8):3324–3331, 2018.
  * [FP09] Georgios E Fainekos and George J Pappas.  Robustness of temporal logic specifications for continuous-time signals.  Theoretical Computer Science, 410(42):4262–4291, 2009.
  * [GKC13] Sicun Gao, Soonho Kong, and Edmund M Clarke.  dreal: An smt solver for nonlinear theories over the reals.  In Automated Deduction–CADE-24: 24th International Conference on Automated Deduction, Lake Placid, NY, USA, June 9-14, 2013. Proceedings 24, pages 208–214. Springer, 2013.
  * [GM20] Maciej Gazda and Mohammad Reza Mousavi.  Logical characterisation of hybrid conformance.  In International Colloquium on Automata, Languages, and Programming (ICALP), 2020.
  * [HWOY22] Timothy G Holloway, John B Williams, Djamila Ouelhadj, and Gong Yang.  Exploring the use of water resource recovery facility instrument data to visualise dynamic resilience to environmental stressors.  Water Research, page 118711, 2022.
  * [IHL23] Roland B Ilyes, Qi Heng Ho, and Morteza Lahijanian.  Stochastic robustness interval for motion planning with signal temporal logic.  In 2023 IEEE International Conference on Robotics and Automation (ICRA), pages 5716–5722. IEEE, 2023.
  * [KP21] Vaibhav Katewa and Fabio Pasqualetti.  Optimal dynamic load-altering attacks against power systems.  In ACC’21, pages 4568–4573. IEEE, 2021.
  * [Pil19] Victoria Yan Pillitteri.  Developing cyber resilient systems: A systems security engineering approach.  National Institute of Standards and Technology, 2:800–160, 2019.
  * [R+76] Walter Rudin et al.  Principles of mathematical analysis, volume 3.  McGraw-hill New York, 1976.
  * [RSHL19] David Rehak, Pavel Senovsky, Martin Hromada, and Tomas Lovecek.  Complex approach to assessing resilience of critical infrastructure elements.  International journal of critical infrastructure protection, 25:125–138, 2019.
  * [Sch98] Alexander Schrijver.  Theory of linear and integer programming.  John Wiley & Sons, 1998.
  * [SDSC17] Moritz Schulze Darup, Rainer Schaich, and Mark Cannon.  How scaling of the disturbance set affects robust positively invariant sets for linear systems.  International Journal of Robust and Nonlinear Control, 27(16):3236–3258, 2017.
  * [SGF20] Adnane Saoud, Antoine Girard, and Laurent Fribourg.  Contract-based design of symbolic controllers for safety in distributed multiperiodic sampled-data systems.  IEEE Transactions on Automatic Control, 66(3):1055–1070, 2020.
  * [SGF21] Adnane Saoud, Antoine Girard, and Laurent Fribourg.  Assume-guarantee contracts for continuous-time systems.  Automatica, 134:109910, 2021.
  * [SJS23] Adnane Saoud, Pushpak Jagtap, and Sadegh Soudjani.  Temporal logic resilience for cyber-physical systems.  In Proceedings of the IEEE Conference on Decision and Control, 2023.
  * [WKK22] Jun Wang, Samarth Kalluraya, and Yiannis Kantaros.  Verified compositions of neural network controllers for temporal logic control objectives.  In 2022 IEEE 61st Conference on Decision and Control (CDC), pages 4004–4009. IEEE, 2022.
  * [WVS20] Ben Wooding, Vahid Vahidinasab, and Sadegh Soudjani.  Formal controller synthesis for frequency regulation utilising electric vehicles.  In Smart Energy Systems and Technologies (SEST), pages 1–6. IEEE, 2020.
  * [ZPV19] Shufang Zhu, Geguang Pu, and Moshe Y Vardi.  First-order vs. second-order encodings for-to-automata translation.  In Theory and Applications of Models of Computation (TAMC’19), pages 684–705. Springer, 2019.
  * [ZTL+17] Shufang Zhu, Lucas M. Tabajara, Jianwen Li, Geguang Pu, and Moshe Y. Vardi.  Symbolic LTLf synthesis.  In Proceedings of the Twenty-Sixth International Joint Conference on Artificial Intelligence, pages 1362–1369, 2017.