2022

Anahita Jamshidnejad*

# Hierarchical Integration of Model Predictive and Fuzzy Logic Control for
Combined Coverage and Target-Oriented Search-and-Rescue via Robots with
Imperfect Sensors

Christopher de Koning<EMAIL_ADDRESS>Control and Operations
Department, Delft University of Technology, Delft, the Netherlands

###### Abstract

Search-and-rescue (SaR) in unknown environments is a crucial task with life-
threatening risks that requires precise, optimal, and fast decisions. Robots
are promising candidates for autonomously performing various SaR tasks in
unknown environments. While humans use their heuristics to effectively deal
with uncertainties of SaR, optimisation of multiple objectives (e.g., mission
time, area covered, number of victims detected) in the presence of physical
and control constraints is a mathematical challenge that requires machine
computations. Thus having both human-inspired and mathematical decision making
capabilities is highly desired for SaR robots, while control approaches that
exhibit both capabilities have been ignored significantly in the literature.
Moreover, coordinating the individual decisions of robots with little
computation cost in large-scale SaR missions is an open challenge. Finally, in
real-life SaR missions due to defects (e.g., in sensors) or due to
environmental factors (e.g., smoke) data perceived by SaR robots may be prone
to uncertainties. We introduce a hierarchical multi-agent control architecture
that exploits non-homogeneous and imperfect perception capabilities of SaR
robots, as well as the computational efficiency and robustness to failure of
decentralised control methods and global performance improvement of
centralised control methods. The integrated structure of the proposed control
framework allows to combine human-inspired and mathematical decision making
methods, via respectively fuzzy logic and model predictive control, in a
coordinated and computationally efficient way. The results of various
computer-based simulations show that while the area coverage of the proposed
approach is comparable to existing heuristic methods that are particularly
developed for coverage-oriented SaR, the efficiency of the introduced approach
in locating the trapped victims is significantly higher. Furthermore, with
comparable computation times, the proposed control approach successfully
avoids potential conflicts that exist in non-cooperative methods. These
results confirm that the proposed multi-agent control system is capable of
combining coverage-oriented and target-oriented SaR in a balanced and
coordinated way.

###### keywords:

Multi-robot search-and-rescue, model predictive control, fuzzy logic control,
imperfect sensors

##  1 Introduction

Search-and-rescue (SaR) robots are expected to take over life-threatening
tasks, especially within initial stages of searching an unknown environment,
in order to reduce the risks for the SaR crew. SaR robots can potentially
reduce the crucial time of finding the trapped victims and allow human
resources to be available for other tasks, e.g., logistics and assisting the
detected victims Casper ; Coburn . Robots can move through areas that are
inaccessible to humans, gather information (e.g., about the location of
victims, explosive materials, and debris) and make maps of the environment.
This way SaR robots contribute to improving the situational awareness for SaR
crews, which is essential for mitigating the mission risks and for saving the
lives of trapped victims Riley_human-robot ; Shimanski2005 ; Chandarana2021 .

SaR approaches can be categorised as target-oriented and coverage-oriented
based on their control objectives. In target-oriented SaR, knowledge about the
target distribution in the environment is initially available (see, e.g.,
Jamshidnejad2018 ; Beck_nonhomocollab ; deAlcantaraAndrade_et_al_2019 ;
sanjuan_fuzzyprioSARS ; Yao_GBNN ). When the SaR environment is unknown,
coverage-oriented approaches are mainly used Galceran_CPPsurvey . Ant colony
algorithms are bio-inspired area coverage methods that are computationally
efficient and easy to implement Koenig_Ants ; Wagner_ants . Machine learning
and neural network methods are also used for area coverage, where robots
progressively learn effective area coverage behaviours Yang_ANNcoverage ;
yang_coopsearch . The main drawback of such methods is their need for training
before they can be implemented. Autonomous learning algorithms, including
generalized model-free reinforcement learning methods, have thus been
developed to address this challenge (see, e.g., Tutsoy2021 ; Tutsoy2017 ). In
these algorithms the system keeps on learning an optimal policy online.
Although very promising, autonomous learning methods may face new challenges
regarding computational burden for real-life implementations for SaR, due to
the large size and varying dynamics of SaR environments that make the learning
procedure more complicated. Moreover, in SaR missions there may be high risks
associated with implementing a solely learning-based algorithm before the
system achieves its optimal performance. More specifically, during the stages
that the algorithm is learning an optimal policy, there are serious risks
regarding losing the trapped victims or delaying their detection, which may
result in their health state becoming critical. Moreover, most coverage-
oriented approaches do not systematically incorporate victim or target
detection in their search behaviour. Arnold et al. in Arnold_nonhomoSARS
present a cooperative, multi-agent SaR system with the objective of both
victim detection and exploration in order to increase situational awareness of
the environment. The SaR agents, however, are steered according to fixed
behaviour sets. This limits the adaptability and thus efficiency of these
robots in highly dynamic SaR environments. Existing SaR control methods are
mainly focused on either coverage or target-oriented SaR. Moreover, MPC, which
is an optimisation-based control method that systematically handles state and
input constraints and can provide robustness to SaR uncertainties, has been
ignored for the crucial task of area coverage in SaR
deAlcantaraAndrade_et_al_2019 . Instead, MPC has mainly been used for
reference tracking in target-oriented SaR in (partially) known environments
(see, e.g., Jamshidnejad2018 ; farrokhsiar_et_al_2013 ; hoy_et_al_2012 ).

In order to speed up mapping the SaR area and to reduce the risk of mission
failure, a fleet of SaR robots may be deployed (see, e.g., cooper_2020 ;
deAlcantaraAndrade_et_al_2019 ; paez_et_al_2021 ). In centralised SAR multi-
agent control, robots are controlled via a centralised system that determines
the mission plans for all these robots (see, e.g., Beck_nonhomocollab ;
deAlcantaraAndrade_et_al_2019 ; sanjuan_fuzzyprioSARS ; Yao_GBNN ). In
decentralised SaR multi-agent control, local (on-board) controllers are
considered for robots (see Koenig_Ants ; Yang_ANNcoverage ; Arnold_nonhomoSARS
; Choi_decentral ; Liu_SARsurvey . Best et al. Best2020  present a cooperative
distributed information gathering approach for SaR robots where based on
learning and heuristics robots visit stationary, pre-known goal regions. While
a task assignment problem is solved in a communication-wise efficient way,
there are no (dynamic) uncertainties involved in the environment of the
robots. Otte et al. Otte2020  address a cooperative task-assignment problem
using a decentralised auction approach. In particular, they investigate the
effect of lossy communication among the agents on the performance of the
multi-agent system with the aim of providing insight into the selection of an
auction algorithm that, despite lossy communication, satisfies the desired
performance criteria of a multi-agent system. A multi-agent search-planning
approach is introduced in Kashino2020  for wilderness SaR with a team of
aerial and ground robots. In their approach, the initial trajectory planning
for the aerial robots is performed offline. After an aerial robot detects a
(possibly moving) target, the robot tracks it until a ground robot intercepts
this target.

While decentralised control approaches are more robust to failure and are
computationally more efficient than centralised approaches, providing reliable
and stable communication among the robots and missing a global vision of the
entire system are challenges of decentralised control approaches
Choi_decentral . Hierarchical architectures can combine the strengths of
centralised and decentralised control approaches (see, e.g.,
Tol_Jamshidnejad_2021 ). Particularly, for multi-robot SaR systems
hierarchical control architectures can provide coordination in the behaviour
of local controllers. However, a limited amount of research on hierarchical
control for SaR robots is available. Examples include Khamis ; Elston_hierarch
; Chandler_hierarch , which are all limited to target-oriented SaR.

Currently, SaR robots need (intensive) supervision and control from human
operators. On the one hand, however, for safety, efficiency, and avoiding
additional challenges regarding online human-robot interaction there is
interest in making these robots autonomous Liu_SARsurvey ; Pellerin_surveySARS
; Casper ; Riley_human-robot . On the other hand, humans use their heuristics
effectively in order to deal with uncertainties of SaR missions. Therefore,
providing SaR robots with human knowledge will improve their performance.
Human knowledge, which is provided as information-based control for SaR robots
in krzysiak_butail_2021 , was shown to improve the performance of SaR robots
in finding the targets. Thus having both capabilities of human-inspired
decision making and mathematical control is highly desired for SaR robots.
However, control approaches that exhibit both capabilities have been ignored
significantly in the literature.

In this paper, we introduce a hierarchical control architecture for multi-
agent control of SaR robots with non-homogeneous, imperfect sensors that
combines mathematical and human-inspired control methods in a computationally
efficient way. The main contributions of this paper include:

  1. 1.

Introducing a novel hierarchical control framework for multi-objective control
and coordination of multi-robot systems and for exploiting their non-
homogeneous sensor imperfections in unknown environments. The resulting
control system benefits both from robustness to failure and computational
efficiency of decentralised control methods and from globally effective
performance of centralised control methods

  2. 2.

Integrating human-inspired and mathematical decision making by formulating
local fuzzy logic controllers, which mimic decision making of human experts,
and a supervisory model predictive control (MPC) system, which provides
mathematical precision in the decisions of the system, systematically handles
state and input constraints, improves the global performance of the multi-
robot system based on its optimal and predictive decision making, and resolves
conflicts of local heuristic controllers

  3. 3.

Implementing the proposed control approaches for combined coverage and target-
oriented SaR via multi-robot systems with imperfect sensors for optimising the
mission time, area coverage, and number of detected victims, and running
extensive experiments via computer-based SaR simulations in order to evaluate
various performance criteria (e.g., computational efficiency, percentage of
area covered, overall certainty level of the map developed for the SaR
environment, number of victims detected) of the proposed SaR control methods
compared to the state-of-the-art methods

Additionally, since the local controllers steer the robots, the multi-robot
control system will not fail due to a failure of the centralised controller.
In that case, as our simulation results indicate, the remaining decentralised
control system can still steer the SaR system safely, although with a degraded
performance. Finally, MPC - which is an optimisation-based control method that
systematically handles state and input constraints and that can provide
robustness to SaR uncertainties - has been ignored in the literature for the
crucial task of area coverage in SaR. Instead, MPC has mainly been used for
reference tracking in target-oriented SaR in (partially) known environments.
Our novel approach and formulation for multi-agent control systems enables MPC
to provide all its strong points for, not only target-oriented, but also
coverage-oriented SaR in unknown areas.

The rest of the paper is structured as it follows. In Section 2 the problem
formulation is detailed. Section 3 discusses the proposed hierarchical mission
planning control approach for SaR robots. Section 4 describes the case study
and experimental setup and presents, analyses, and discusses the results.
Section 5 concludes the paper and provides suggestions for future research.

##  2 Problem Formulation

In this section, we explain and formulate the details of the mission planning
problem of SaR robots with non-homogeneous imperfect sensors. In particular,
we discuss the modelling of the SaR environment and victims and the
uncertainties involved, as well as the formulation of perception capabilities
of SaR robots.

###  2.1 SaR Environment

Figure 1: Schematic view of a SaR environment.

The SaR environment $E$ is modelled by a bounded, discretised, 2-dimensional
cellular area of $L_{x}\times L_{y}$ cells (see Figure 1.). Each cell in the
SaR environment corresponds to the coordinates $(x,y)$ of its centre and may
be empty or occupied by a static obstacle (i.e., wall, pillar, rubble), or by
a victim and/or a SaR robot. A cell can embed a single victim at a time.
Moreover, obstacles make a cell inaccessible for SaR robots and for victims.
The following uncertainties for SaR robots have been considered:

  * •

External uncertainty regarding the SaR environment, i.e., the total number of
victims and obstacles and their positions are unknown

  * •

External, random uncertainty about the pattern of movement of victims

  * •

Internal (i.e., structural) and external (i.e., proximal) uncertainties
regarding the perceived data

An $L_{x}\times L_{y}$ matrix $\mathcal{W}(\kappa)$, called the _occupancy
map_ , is used to record cells that are occupied by static obstacles after
being detected by a SaR robot. Furthermore, whenever a victim is detected by
the sensor of SaR robot $i$, the robot stores the location, perceived health
state, and time of detection of the victim in a local matrix (specific to the
robot) called the _victim map_ $\mathcal{V}_{i}(\kappa)$ of SaR robot $i$.
This map is used by the controller of the robot to make the current control
decision. However, for the sake of efficiency for the on-board computations
and the memory storage, in their local victim maps the robots keep track only
of those victims who have been selected as a target by the controller and have
been visited by the robot (i.e., the robot has been in the same cell as the
victim), as well as of those victims who are currently within the perception
field of the robot’s sensor. Thus SaR robots do not record any memory of
victims who have previously been detected by the robot, but have not been
selected as target for them.

The scan certainty $c(x,y,\kappa)$ of cell $(x,y)$ for time step $\kappa$ is a
value within $[0,1]$ that specifies the certainty level regarding the
information available about cell $(x,y)$. Moreover, each cell $(x,y)$
corresponds to a proximal uncertainty at time step $\kappa$ that is a function
of the Euclidean distance of the cell from all SaR robots that scan the cell.
The scan certainty of a cell depends on whether or not the cell has been
scanned by any SaR robots and if so, how accurate the perceived data is, i.e.,
the scan certainty depends on the proximal uncertainty. This relationship is
explained in detail in Section 2.2. Initially the scan certainty of all cells
within the SaR environment is zero. The scan certainty for all cells is
included in an $L_{x}\times L_{y}$ matrix $\mathcal{C}(\kappa)$, called the
_scan certainty map_ , which will be updated in time.

###  2.2 SaR Robots

We consider a multi-robot SaR system composed of $N$ agents $\textrm{a}_{i}$
($i=1,\ldots,N$) that, per simulation time step, may move to one of the $8$
neighbouring cells (see Figure 2(a)). These robots are equipped with optical
cameras and sensors that localise the victims and that assess their health
state (e.g., acoustic and heat sensors Casper ; Ganesan  or sensors that
detect WiFi-enabled devices Wang_wifiSARS ). The perception field
$E_{i}(\kappa)$ of SaR robot $i$ for time step $\kappa$ includes all cells of
the SaR environment that fall within a circle of radius $r_{\mathrm{p},i}$,
centred at the position of the robot at time step $\kappa$, where
$E_{i}(\kappa)\subseteq E$ (see Figure 2(b)).

The data perceived by SaR robots may in general be imperfect, i.e., scanning a
cell does not necessarily yield full certainty about the information within
the cell. Two sources of uncertainty regarding the perceived data are
considered: (1) _Structural imperfection_ , which corresponds to a fixed
perceptual uncertainty reduction rate $\eta_{i}\in(0,1]$ per SaR robot
$\textrm{a}_{i}$. More specifically, every time SaR robot $\textrm{a}_{i}$
scans a cell, the uncertainty regarding the information of the cell is reduced
by rate $\eta_{i}$. Thus when $\eta_{i}=1$ there is no structural imperfection
(i.e., the perceived information corresponds to $100\%$ certainty). Moreover,
we do not consider sensors that are completely out of function (i.e.,
$\eta_{i}=0$) due to structural imperfection. (2) _Proximal uncertainty_ ,
which implies that while all cells within the perception field $E_{i}$ of SaR
robot $\textrm{a}_{i}$ are scanned, the degree of increase in the scan
certainty of these cells decreases according to their distance from the
sensor.

The structural imperfection of sensors and the proximal uncertainty of the
cells together will result in an uncertainty dynamic ratio
$\sigma(x,y,\kappa)$ corresponding to every cell $(x,y)$ per time step
$\kappa$. We have:

$\displaystyle z(x,y,\kappa+1)=\sigma(x,y,\kappa)z(x,y,\kappa)$ (1)

with $z(x,y,\kappa)$ the scan uncertainty (i.e., $1-c(x,y,\kappa)$) assigned
to cell $(x,y)$ at time step $\kappa$, and:

$\displaystyle\sigma(x,y,\kappa)$
$\displaystyle=\prod_{i=1}^{N}\sigma_{i}(x,y,\kappa)$ (2)
$\displaystyle\sigma_{i}(x,y,\kappa)$
$\displaystyle=1-(1-\eta_{i})e^{-r_{i}(x,y,\kappa)}\cdot$ (3)
$\displaystyle\frac{\displaystyle
1-\operatorname{sign}\Big{(}r_{i}(x,y,\kappa)-r_{\mathrm{p},i}\Big{)}}{\displaystyle
2}$

where $r_{i}(x,y,\kappa)$ is the Euclidean distance of SaR robot
$\textrm{a}_{i}$ to cell $(x,y)$ at time step $\kappa$,
$\operatorname{sign}(\cdot)$ represents the sign function, and
$\sigma_{i}(x,y,\kappa)$ is the share of the uncertainty dynamic ratio of cell
$(x,y)$ at time step $\kappa$ that is provided by the sensor of SaR robot
$\textrm{a}_{i}$. The updated scan certainty for cell $(x,y)$ is given by:

$\displaystyle c(x,y,\kappa+1)=1-z(x,y,\kappa+1)$ (4)

Based on (3), the effect of the proximity on the uncertainty dynamic ratio of
cells is modelled by an exponential function. More specifically, when
$r_{i}(x,y,\kappa)=0$, i.e., for the cell where SaR robot $\textrm{a}_{i}$ is
currently located at, the uncertainty dynamic ratio corresponding to SaR robot
$\textrm{a}_{i}$ is $\eta_{i}$ (i.e., the maximum possible improvement in the
scan certainty of the cell that can be provided by the sensor of SaR robot
$\textrm{a}_{i}$). This uncertainty dynamic ratio varies exponentially until
for $r_{i}(x,y,\kappa)\geq r_{\mathrm{p},i}$, it becomes unity (i.e., the scan
certainty of the cell at the current time step does not improve as a result of
a contribution of the sensor of SaR robot $\textrm{a}_{i}$).

SaR robots may differ from each other in two properties regarding their
sensors: (1) The sensors of SaR robots may have different perception radii
$r_{\mathrm{p},i}$ for $i=1,\ldots,N$. (2) The accuracy of these sensors, and
thus their perceptual uncertainty reduction rate $\eta_{i}$ may be different.

Figure 2: (a) Movement possibilities of a SaR robot. (b) Perception field of a
SaR robot.

###  2.3 Victim Modelling

The number, location, and health state of the victims are initially unknown
for the SaR robots. The victims follow a random pattern of movement, i.e.,
victim $v$ with position
$(x^{\textrm{v}}_{v}(\kappa),y^{\textrm{v}}_{v}(\kappa))$ at time step
$\kappa$ may remain in its current cell with probability $p^{\textrm{s}}$ or
may move to one of its (unblocked) neighbouring cells with a total probability
$1-p^{\textrm{s}}$, which results in an equal probability of
$p^{\textrm{m}}_{v}(\kappa)=(1-p^{\textrm{s}})/n^{\textrm{free}}(x^{\textrm{v}}_{v}(\kappa),y^{\textrm{v}}_{v}(\kappa),\kappa)$
to move to each of the (unblocked) neighbouring cells, where
$n^{\textrm{free}}(x,y,\kappa)$ is the number of free neighbouring cells for
cell $(x,y)$ at time step $\kappa$. Note that for the sake of simplicity we
consider the probability $p^{\textrm{s}}$ to be constant in time, space, and
for all victims. In case a more detailed model is to be used, this probability
may vary with time, and per cell and victim.

Moreover, each victim holds a certain health state, $h_{v}(\kappa)$, which
varies within $[0,100]$ and implies how healthy or injured the victim is at
time step $\kappa$. Whenever a victim is detected by a SaR robot, their
initial health state is registered. Over time, the health state of each victim
may decrease with the rate $\Delta h_{v}(\kappa)$ given by:

$\displaystyle{\Delta}h_{v}(\kappa)=\bigg{\\{}\begin{array}[]{cc}-\alpha&h^{\textrm{crit}}\leq
h_{v}(\kappa)\leq 100\\\ \beta h_{v}(\kappa)-\gamma&0\leq h_{v}(\kappa)\leq
h^{\textrm{crit}}\end{array}$ (7)

with $\alpha,\beta,\gamma>0$, $h^{\textrm{crit}}$ the critical health state,
and $\gamma\geq\beta h^{\textrm{crit}}$. Based on (7), a victim has a
uniformly deteriorating health state whenever their health state is not less
than $h^{\textrm{crit}}$ (i.e., health state is stable), while the rate of
deterioration of the health state becomes linear as soon as the health state
is below $h^{\textrm{crit}}$. The updated health state is given by:

$\displaystyle
h_{v}(\kappa+1)=\max\left\\{h_{v}(\kappa)+{\Delta}h_{v}(\kappa),0\right\\}$
(8)

A SaR robot detects a victim whenever they are both in the same cell. Without
considering the technical details regarding data analysis, sensor fusion, or
soft sensing in this paper, we assume that the robot detects the victim (see,
e.g., Dousai2022 ; Llasag2019 ) and assesses the health state of the victim,
e.g., using a combination of WiFi, optical, thermal, or acoustic sensors and
using image processing algorithms or via direct feedback received from the
victims when possible (see, e.g., Pinheiro2022 ).

##  3 Hierarchical Control System

Figure 3: Hierarchical architecture of the proposed cooperative hierarchical
mission planning controller.

Next we explain the proposed hierarchical control system that steers the
search behaviour of SaR robots. The control architecture includes two levels
(see Figure 3): The lower level of control is composed of decentralised
controllers that steer the local search behaviour of each SaR robot, while the
higher control level includes a centralised supervisory controller that
coordinates the behaviour of the decentralised controllers, such that search
conflicts among SaR robots are resolved. SaR robots only communicate with the
supervisory control level, without sharing any information among themselves.
The proposed control architecture thus combines the strengths of centralised
and decentralised control approaches.

###  3.1 Local Fuzzy Logic Controllers

At the local level, a SaR robot first processes the data that is captured via
its sensors and then constructs a local priority map for its perception field.
This map includes quantities corresponding to the importance of visiting the
cells for the SaR robot. Next the local controller of the robot determines a
path that yields the highest local gain according to a quantity called the
path grade. Since the main objective of the SaR mission is to optimise the
area coverage and time efficiency of detecting the victims, the following two
main criteria are considered in grading a path:

  1. 1.

Time reduction: Each SaR robot should reach its targets in the least possible
time in order to contribute to minimising the overall mission time.

  2. 2.

Exploration increase: Each SaR robot should scan as many (unexplored) cells as
possible along its path in order to contribute to maximising the overall area
coverage.

These two criteria may possess a conflict, since for the first criterion the
search behaviour should be target-oriented (in other words the robot should
find the shortest possible path that leads it to the target as quickly as
possible), while for the second criterion the search behaviour is coverage-
oriented (in other words the robot should visit more cells before reaching its
target). Therefore, the local controllers are developed such that a balanced
trade-off between these criteria is provided.

####  3.1.1 Search Priority Assignment

The local controllers of SaR robots should first assign priorities to
potential paths, specifying the urgency to scan every cell within the
environment. For every cell in the perception field $E_{i}$ of SaR robot
$\textrm{a}_{i}$, a priority score is determined using rule-based fuzzy logic
control (FLC) methods. The main motivation for using FLC is its computational
efficiency, which is essential for local controllers due to the limited
computational power available on board for SaR robots, and the capability of
FLC in mimicking human’s logic for decision making, which allows to
incorporate human expert knowledge in local controllers. Thus local
controllers effectively mimic the reasoning of human experts without their
direct supervision.

Fuzzy rules with the following formulation are used by local controllers of
SaR robots:

$\displaystyle\mathcal{R}_{m}:$ $\displaystyle\text{{ If
}}e^{\textrm{v}}(x,y,\kappa)\text{ is }A_{m,1}\text{{ and
}}h_{v}(\kappa)\text{ is }A_{m,2}$ (9) $\displaystyle\text{{ and
}}c(x,y,\kappa)\text{ is }A_{m,3}\text{{ then }}\rho(x,y,\kappa)\text{ is
}B_{m}$

with $e^{\textrm{v}}(x,y,\kappa)$ the probability of existing a victim in cell
$(x,y)\in E_{i}$ at time step $\kappa$, $m=1,\ldots,M$ with $M$ the number of
rules, $A_{m,1}$, $A_{m,2}$, $A_{m,3}$, and $B_{m}$ fuzzy sets that adopt a
linguistic term.

The fuzzy inference system corresponding to the rules given by (9) receives
$3$ inputs (i.e., the probability of existence of a victim in a cell, the
health state of the potential victim, and the most recent scan certainty value
of a cell) and assigns a search priority $\rho(x,y,\kappa)$ to cell $(x,y)$
for time step $\kappa$. Note that every SaR robot $\textrm{a}_{i}$ has access
to its local knowledge stored in the local scan certainty map
$\mathcal{C}_{i}(\kappa)$ and local victim map $\mathcal{V}_{i}(\kappa)$. The
probability $e^{\textrm{v}}(x,y,\kappa)$ of existence of a victim in cell
$(x,y)$ estimated by SaR robot $\textrm{a}_{i}$ depends on the robot’s sensor,
and adopts either a very small value when the sensor receives no signal that
implies a victim exists in the cell (for a sensor with $\eta_{i}=1$ this small
value may be $0$, while for a sensor with $\eta_{i}\in(0,1)$ a small positive
value may be considered), or a percentage determined according to the
structural imperfection and proximal uncertainty explained in Section 2.2.
Based on (9), cells within the perception field of the robot where it is
likely to find a victim with a worse health state and has not yet
(extensively) been scanned receive a higher priority. Inaccessible cells
within the occupancy map $\mathcal{W}(\kappa)$ receive a null priority.

####  3.1.2 Path Planning

After prioritising the cells, each local controller determines potential paths
for the corresponding SaR robot. In order to optimise the time, shortest paths
are favourable, while for optimising the area coverage, paths that visit more
cells with higher priorities are preferred. Thus the local controller applies
an A* search approach Hart_ASTAR  based on Yen’s algorithm Yen_kshortest  to
determine a certain number of shortest paths that end at every cell within the
perception field of the robot. Afterwards these paths are graded based on
their travel time and degree of exploration to specify how favourable they are
for the SaR mission at the current time step. The travel time is computed
based on the path length and the robot’s speed. We suppose that a SaR robot
moves one cell per time step, thus the travel time corresponds to the path
length only. The degree of exploration of every potential path
$P_{i}(\kappa)\subseteq E_{i}$ for SaR robot $\textrm{a}_{i}$ at time step
$\kappa$ is computed via:

$\displaystyle\epsilon\left(P_{i}(\kappa)\right)=\sum_{k=\kappa}^{\kappa+\ell\left(P_{i}(\kappa)\right)-1}\lambda^{k}\rho(x^{\textrm{a}}_{i}(k),y^{\textrm{a}}_{i}(k),k)$
(10)

where the path is defined by:

$\displaystyle P_{i}$
$\displaystyle(\kappa)=\Big{\\{}\left(x^{\textrm{a}}_{i}(\kappa),y^{\textrm{a}}_{i}(\kappa)\right),\ldots,$
(11)
$\displaystyle\left(x^{\textrm{a}}_{i}\left(\kappa+\ell(P_{i}(\kappa))-1\right),y^{\textrm{a}}_{i}\left(\kappa+\ell(P_{i}(\kappa))-1\right)\right)\Big{\\}}$

with $\ell(P_{i}(\kappa))$ the path length, $\lambda\in[0,1]$ the discount
factor, and $\rho(x^{\textrm{a}}_{i}(k),y^{\textrm{a}}_{i}(k),k)$ the priority
value of cell $(x^{\textrm{a}}_{i}(k),y^{\textrm{a}}_{i}(k))$, which the robot
should visit at time step $k=\kappa,\ldots,\kappa+\ell(P_{i}(\kappa))-1$ when
it follows path $P_{i}(\kappa)$. Note that since the priority values for the
cells corresponding to time steps $k>\kappa$ are based on predicted estimates,
considering a discount factor can reduce the potential influence of errors in
the predictions. Finally, the grade of path $P_{i}(\kappa)$ is computed by
(with $c_{1},c_{2}>0$ constant values):

$g\left(P_{i}(\kappa)\right)=-c_{1}\ell\left(P_{i}(\kappa)\right)+c_{2}\epsilon\left(P_{i}(\kappa)\right)$
(12)

######  Remark 1.

Since the paths that will be generated by the control system of the SaR robots
are rectilinear, for practical implementations and to make the path easy to
execute for real robots, we propose smoothening the paths before
implementation (see, e.g., Jamshidnejad_DeSchutter2018  for equations that can
be used to smoothen such paths).

###  3.2 Supervisory MPC Controller

At the supervisory level, a centralised MPC-based controller is used that
receives the local information corresponding to each SaR robot and merges them
to build up global maps of the current perception fields of the robots. Note
that while robots erase the non-target victims from their local victim maps
(see Section 2.1 for details), the global victim map keeps track of all
locally perceived information. This is practically possible because the global
maps are recorded on a remote computer station that is not restricted by
computational and memory limits.

The supervisory controller is called whenever a search conflict is identified,
i.e., the cardinality of the intersection of the perception fields of two SaR
robots exceeds a certain threshold: $\textrm{card}\left(E_{i}\cap
E_{j}\right)>\tau_{\textrm{int}}$. A model of the environment including the
most updated cognitive maps is used as the prediction model of the supervisory
controller, which determines globally optimal (within the controller’s
prediction time window) paths for SaR robots. This optimality is defined as a
trade-off among various objectives including the mission time, the area
coverage, and the chances of visiting more trapped victims with a more crucial
health state. Despite providing globally optimal solutions, the MPC controller
is computationally demanding due to the size of the centralised optimisation
problem and the non-linearities involved in the problem. Therefore, we provide
the supervisory controller with the paths that are determined by the local
controllers as a warm start for the MPC optimisation problem to converge
faster to an optimal solution. Taking into account the objectives of the SaR
mission, the objective function to be maximised by the supervisory controller
at time step $\kappa$ is given by:

$\displaystyle J\left(\mathbb{P}(\kappa)\right)=$ $\displaystyle
w_{1}\sum_{i=1}^{N}g\left(P_{i}(\kappa)\right)+w_{2}\sum_{(x,y)\in
E}c\left(x,y,N^{\textrm{p}}(\kappa)\right)$ (13)

with $\mathbb{P}(\kappa)$ (the optimisation variable) the set of paths for all
the $N$ SaR robots and $w_{1}$ and $w_{2}$ constant weights. The objective
function given by (3.2) is a weighted sum of two terms: (i) the overall grade
of all paths (estimated by (12)) and (ii) the total predicted scan certainty
of the SaR environment at the end of the current prediction horizon
$N^{\textrm{p}}(\kappa)$, which is given by
$N^{\textrm{p}}(\kappa)=\max_{i=1,\ldots,N}{\ell(P_{i}(\kappa))}$. Thus the
second term steers the fleet of the SaR robots to spread out over the
environment. In other words, the supervisory controller provides a balanced
trade-off between locally preferred paths per robot and globally optimal paths
from the point of area coverage. The supervisory controller does this using a
global scan certainty map $\mathcal{C}(\kappa)$ of the environment and a
global victim map $\mathcal{V}(\kappa)$, which are built by merging the local
maps of SaR robots.

The supervisory control optimisation problem for time step $\kappa$ is given
by (where the prediction window is
$\\{\kappa,\ldots,N^{\textrm{p}}(\kappa)-1\\}$):

$\displaystyle\max_{\mathbb{P}(\kappa)}J(\mathbb{P}(\kappa))$ such that:
$\displaystyle\mathbb{P}(\kappa)=\left\\{P_{1}(\kappa),\ldots
P_{N}(\kappa)\right\\}$ (14a) For all the paths, (11) holds, with (14b)
$\displaystyle\hskip
64.58313pt\left(x^{\textrm{a}}_{i}(\kappa),y^{\textrm{a}}_{i}(\kappa)\right)\in
E\setminus\mathcal{W}(\kappa)$
$\displaystyle\left(x^{\textrm{v}}_{v}(\kappa),y^{\textrm{v}}_{v}(\kappa)\right)\notin
P_{i}(k)\cap P_{j}(k),\textrm{ where}$ (14c) $\displaystyle\hskip
30.1388pti,j=1,\dots N,i\neq j,v=1,\ldots,N^{\textrm{v}}(\kappa)$
$\displaystyle\left(x^{\textrm{a}}_{i}(k),y^{\textrm{a}}_{i}(k)\right)=\left(x^{*}_{i}(\kappa),y^{*}_{i}(\kappa)\right)\qquad
i=1,\dots N$ (14d)

Constraints (14a) and (14b) define the optimisation variable and state that
the paths should be feasible. Constraint (14c) restricts multiple SaR robots
to visit the same victim, where $N^{\textrm{v}}(\kappa)$ is the number of
victims detected until simulation time step $\kappa$. This constraint improves
the victim search efficiency and area coverage of the robots. To reduce the
conservativeness of the problem and avoid infeasibility, constraint (14c) may
be defined as a chance (instead of a hard) constraint. Finally, constraint
(14d) allows the starting point of the paths to be the most recent measured
coordinates $(x^{*}_{i}(\kappa),y^{*}_{i}(\kappa))$ of the corresponding SaR
robot.

######  Remark 2.

Since the objective function of the supervisory MPC-based controller is
defined in (3.2) as a weighted sum of the multiple control objectives, these
objective terms will be normalised when implementing the optimisation problem.

##  4 Case Study

Next we discuss the results of computer-based simulations that are
systematically designed to evaluate the performance of the proposed
hierarchical control approach in comparison with state-of-the-art approaches
for SaR. The simulations are implemented via MATLAB R2019b on a PC with Intel
Core i7 Processor with $2.20$ GHz frequency. Whenever an optimization problem
should be solved to determine the paths of the SaR robots, the path planning
problem is solved using pattern search as optimization method, since this
algorithm showed to be faster than other alternative approaches. In order to
make sure that the resulting paths meet the requirements of a discrete
cellular environment for the numerical simulations (i.e., the way points
defining the path have to be located at the centre of a cell) the continuous
coordinates for the way points determined by pattern search are projected to
the centre of the cells using the round function. For the parameters of the
algorithms, we did a manual tuning with respect to the default settings.

###  4.1 Simulation Setup

We consider the following four search approaches that are common for SaR and
compare their performance, in terms of victim detection, area coverage, and
computational efficiency, with the proposed hierarchical control approach,
which we call cooperative controller due to the supervisory MPC level. Selfish
controller: A control system composed of the local controllers described in
Section 3.1, where the main difference with the cooperative controller is the
lack of a supervisory controller. These controllers make decisions that fit
their own circumstances only. Pure MPC controller: An optimisation-based
search approach with the MPC structure of the supervisory controller described
in Section 3.2, where the main difference with the cooperative controller is
the lack of warm starting with trajectories that are proposed by the local
controllers. Instead, as it is common in the implementation of MPC, the MPC
controller receives the shifted solution of the previous time step as a warm
start (see Diehl2005  for more details). Note that, in order to account for
non-convexity of the problem, we have also run the simulations for pure MPC
with multiple starting points within the given time budget. However, the
results for pure MPC with warm start were still better. Thus only the results
regarding pure MPC with warm start have been presented in the paper.

ACS controller: A heuristic ant-colony-based search approach based on
Koenig_Ants , where the global scan certainty map $\mathcal{C}(\kappa)$ is
used for pheromone map for the ant colony system.

Exhaustive controller: A random search strategy for SaR robots, where such
search strategies are commonly used as reference bases for the other search
methods.

Table 1: Modelling and control parameters Parameter | Value |
---|---|---
$L_{x}$ | 40 |
$L_{y}$ | 25 |
Number of victims | 25 |
$p^{\textrm{s}}$ | 0.6 |
$\alpha$ | 0.25 |
$\beta$ | 1/60 |
$\gamma$ | 1 |
$h^{\textrm{crit}}$ | 30 |
$\lambda$ | 0.6 |
$c_{1}$ | 2.0 |
$c_{2}$ | 5.0 |
$\tau_{\textrm{int}}$ | 30 |
$w_{1}$ | 1.0 |
$w_{2}$ | 0.05 |
|  |
Table 2: Rule base of the fuzzy inference system (with the inputs given by
$e^{\textrm{v}}$, i.e., the probability of existence of a victim, $h_{v}$,
i.e., the health state of the potential victim, and $c$, i.e., the scan
certainty of the cell to be visited, and the output is $\rho$, i.e., the
search priority of the cell)

$\mathcal{R}_{m}$ | $e^{\textrm{v}}$ | $h_{v}$ | $c$ | $\rho$
---|---|---|---|---
1 | Low | Stable | Known | Very Low
2 | Low | Medium | Known | Very Low
3 | Low | Stable | Partial | Very Low
4 | Low | Medium | Partial | Low
5 | Low | Critical | Known | Low
6 | Medium | Medium | Partial | Low
7 | Medium | Critical | Known | Low
8 | Low | Stable | Unknown | Low
9 | Medium | Stable | Known | Medium
10 | Medium | Medium | Known | Medium
11 | Medium | Stable | Partial | Medium
12 | High | Stable | Partial | Medium
13 | High | Medium | Known | Medium
14 | High | Critical | Known | Medium
15 | Low | Critical | Partial | Medium
16 | Low | Medium | Unknown | Medium
17 | High | Stable | Known | High
18 | Medium | Stable | Unknown | High
19 | Medium | Medium | Unknown | High
20 | High | Stable | Unknown | High
21 | Low | Critical | Unknown | High
22 | Medium | Critical | Partial | Very High
23 | High | Medium | Partial | Very High
24 | Medium | Critical | Unknown | Very High
25 | High | Medium | Unknown | Very High
26 | High | Critical | Partial | Very High
27 | High | Critical | Unknown | Very High

Table 3: Parameters of SaR robots

General case | $r_{\mathrm{p},i}$ | $\eta_{i}$ | $(x^{\textrm{a}}_{i,0},y^{\textrm{a}}_{i,0})$ |  | Case 3 | $r_{\mathrm{p},i}$ | $\eta_{i}$ | $(x^{\textrm{a}}_{i,0},y^{\textrm{a}}_{i,0})$
---|---|---|---|---|---|---|---|---
$i=1$ | 6 | 0.1 | (1,16) |  | $i=1$ | 7 | 0.1 | (9,10)
$i=2$ | 4 | 0.3 | (13,25) |  | $i=2$ | 3 | 0.3 | (5,8)
|  |  |  |  |  |  |  |
Case 1 | $r_{\mathrm{p},i}$ | $\eta_{i}$ | $(x^{\textrm{a}}_{i,0},y^{\textrm{a}}_{i,0})$ |  | Case 4 | $r_{\mathrm{p},i}$ | $\eta_{i}$ | $(x^{\textrm{a}}_{i,0},y^{\textrm{a}}_{i,0})$
$i=1$ | 6 | 0.1 | (10,8) |  | $i=1$ | 4 | 0.1 | (8,7)
$i=2$ | 4 | 0.3 | (10,6) |  | $i=2$ | 4 | 0.3 | (8,9)
|  |  |  |  |  |  |  |
Case 2 | $r_{\mathrm{p},i}$ | $\eta_{i}$ | $(x^{\textrm{a}}_{i,0},y^{\textrm{a}}_{i,0})$ |  |  | Case 5
---
$r_{\mathrm{p},i}$ | $\eta_{i}$ | $(x^{\textrm{a}}_{i,0},y^{\textrm{a}}_{i,0})$
$i=1$ | 6 | 0.1 | (6,10) |  | $i=1$ | 7 | 0.1 | (9,17)
$i=2$ | 4 | 0.3 | (6,8) |  | $i=2$ | 3 | 0.3 | (5,15)

Figure 4: Trapezoidal fuzzy membership functions defined for the inputs and
outputs of the Mamdani rule bases.

######  Remark 3.

Note that, ideally, a centralised MPC controller can provide the desired
performance for a system by providing a globally optimal solution. This
requires to provide enough computational resources and time for the
centralised MPC controller. However, a main challenge that needs to be
addressed for SaR problems is to provide a balanced trade-off between
performance and computation time, such that the developed control system can
meet the real-time requirements of a SaR robotic team. Therefore, we are
interested in assessing how well different control approaches can steer the
behaviour and performance of the SaR system when they are constrained by the
computation time. Thus for both the pure MPC controller and the supervisory
MPC controller we have considered a limited time budget, which may in some
cases imply a degradation of the performance to meet the given computation
time.

A set of $20$ simulation scenarios, each lasting $300$ simulation time steps,
with a seeded random placement of victims and obstacles in an environment of a
fixed size is considered. The parameters required for these simulations to
estimate the movement of the victims and to compute (7) and (10) are given in
Table 1. The coefficients/weights in (12) and (3.2) are also given in Table 1,
where the corresponding values are tuned manually via extensive experiments.
The terms that describe the sets $A_{m,1}$, $A_{m,2}$, $A_{m,3}$, and $B_{m}$
in (9) should for real-life scenarios be deduced from real expert knowledge.
For the numerical simulations designed in our case studies, we have defined
the corresponding rule base based on intuition. More specifically, the sets
$A_{m,1}$, $A_{m,2}$, and $A_{m,3}$ are verbally described by, respectively,
“Low, Medium, High”, “Critical, Medium, Stable”, and “Unknown, Partially
Known, Known”. This selection allows us to build up a Mamdani rule base
composed of $27$ rules. For the output set $B_{m}$ we select one of the
following terms, “Very Low”, “Low”, “Medium”, “High”, “Very High”, based on
intuition and suited for the given realisations of the input fuzzy sets. The
resulting Mamdani rule base is represented in Table 2. The corresponding
membership functions used in (9) are shown in Figure 4, where trapezoidal
functions have been selected, since they have proven to result in good quality
control systems in various real-life applications (see FMF  for details).

For the case study, we consider $2$ SaR robots with different sensory
perception radii and perceptual uncertainty reduction rates. The robots start
at fixed coordinates without any initial information about the SaR
environment, thus all maps are initialised to zero/null. The parameters used
for the SaR robots are shown in Table 3. These parameters have been selected
such that the influence of non-homogeneity and imperfection of the sensors can
properly be incorporated into the numerical simulations.

In order to evaluate various search approaches in terms of area coverage,
victim search efficiency, and computational efficiency, the following
performance metrics are considered. The area coverage is assessed via two
performance metrics: (1) Total scan certainty of the environment as a function
of simulation time steps, i.e.:

$S(\kappa)=\sum_{\forall(x,y)\in E}c(x,y,\kappa)$ (15)

(2) Rise time for the total scan certainty (of a particular percentage). These
performance metrics quantify the absolute area coverage, as well as the speed,
thus efficiency, of each search approach. The victim search efficiency is
evaluated via three metrics: (1) Number of (live and deceased) victims
detected per simulation time step. (2) Simulation time step for which each
victim is detected. (3) Health state of each victim at the time of detection.
Finally, the average time for making control decisions per simulation time
step is used to report the computational effort of each control approach.

###  4.2 Results & Discussions

Figure 5: Total scan certainty corresponding to different controllers as a
percentage of the maximum scan certainty that can be obtained for the SaR
environment. Figure 6: Victim detection efficiency for each control approach
evaluated by (a) the number of victims found, (b) the number of victims
deceased, (c) time of detection, and (d) health state of victims when
detected. Table 4: Number of simulation time steps needed to reach certain
degrees of scan certainty.

Degree of scan certainty | $50\%$ | $70\%$ | $80\%$ | $85\%$ | $90\%$
---|---|---|---|---|---
Cooperative controller | 97 | 157 | 208 | 291 | -
Selfish controller | 104 | 167 | 238 | - | -
ACS controller | 97 | 159 | 218 | 249 | 292
Pure MPC controller | 112 | 229 | - | - | -
Exhaustive controller | - | - | - | - | -
|  |  |  |  |

Table 5: Average computation time for decision making per simulation time step
for different controllers.

Cooperative | Selfish | ACS | Pure MPC | Exhaustive
---|---|---|---|---
controller | controller | controller | controller | controller
4.5 [s] | 3.2 [s] |  0 [s] | 8.6 [s] |  0 [s]
|  |  |  |

Figure 7: Number of simulation time steps when conflict is detected for
cooperative and selfish controllers.

Figure 5, Table 5, Figure 6, and Table 5 represent the results of the
simulations including, respectively, the total scan certainty, the number of
simulation time steps required per control approach in order to reach a total
scan certainty of $50\%,70\%,80\%,85\%$, and $90\%$, the victim detection
efficiency, and the average decision making time. Moreover, Figure 7 shows the
number of simulation time steps when search conflicts is registered for the
cooperative and selfish controllers. For the cooperative controller, the
number of registered conflicts is the number of times the supervisory MPC
controller has been activated.

Figure 5 shows that the selfish and cooperative controllers achieve a
comparable total scan certainty of above $80\%$ by the end of the simulations
(to be more accurate, with a slightly higher percentage (around $85\%$) for
the cooperative controller). However, Table 5 shows that the cooperative
controller is significantly faster in reaching particular levels of scan
certainty in later stages of the simulation (e.g., needs $14.4\%$ less time to
reach an $80\%$ total scan certainty). The ACS controller (see Figure 5)
reaches an overall scan certainty that is $6.71\%$ larger than that of the
cooperative controller, with a comparable rise time in earlier stages of
simulation. This is because the objective function of the ACS controller
solely considers the scan certainty map to determine the most favourable next
step for each SaR robot per time step. Thus the ACS controller has a single,
non-conflicting objective as opposed to the selfish and cooperative
controllers. Based on Figure 5 and Table 5, the pure MPC controller performs
much better than the exhaustive controllers, but worse than the cooperative,
selfish, and ACS controllers. This is mainly because the pure MPC approach is
prone to falling within a sequence of local optima since it relies on the
solution of the previous time steps as warm start. This not only highlights
the importance of using more systematic (i.e., in line with the global
objectives of the SaR mission) warm starts for the MPC controller, but also
stresses that the exploratory nature of the selfish and cooperative
controllers plays an important role in avoiding such issues. Finally, the
exhaustive controller shows the worst performance regarding the area coverage,
due to its lack of systematic search objective.

Figures 6(a) and 6(b) show that the cooperative and selfish controllers
achieve a similar number of $25$ victims found and using these controllers
correspond to the least number of victims deceased. However, considering all
the simulation runs, the cooperative controller has lower variances, i.e.,
$0.134$ and $0.345$, for respectively the number of victims found and deceased
than those (i.e., $3.10$ and $2.22$) of the selfish controller. This indicates
a more consistently satisfactory performance for victim search using the
cooperative controller. Both the selfish and cooperative controllers
outperform the pure MPC, ACS, and exhaustive controllers in terms of victim
detection, with the ACS and the exhaustive controller showing the worst
performance (see Figures 6(a) and 6(b)). This is because systematic victim
detection is not an objective for these controllers. More specifically, with
the ACS and exhaustive controllers SaR robots may detect the victims randomly.
The fact that the ACS controller detects more victims than the exhaustive one
is an indirect influence of its higher area coverage (see Figure 5). While the
pure MPC method outperforms both the ACS and exhaustive controllers, compared
to the cooperative and selfish controllers less victims are detected and more
victims are deceased. This is due to the lower area coverage by the pure MPC
controller, which has a negative impact on the victim detection efficiency.
Based on Figures 6(c) and 6(d) the cooperative and selfish controllers perform
equally well considering the detection time and health state of victims. While
the MPC controller detects less victims than the cooperative and selfish
controllers, the detection time and health state of the victims found is at a
comparable level (i.e., it is slightly worse) as those of the cooperative and
selfish controllers.

Finally, the computation time per decision (see Table 5) for the cooperative
controller (i.e., $4.5$ s) is almost half of the computation time when only
MPC is used (i.e., pure MPC controller) to steer the system. Moreover,
compared to the average decision making time of the selfish controller (i.e.,
$3.5$ s), and considering the significantly better performance of the
cooperative controller, this controller is the best choice among all the given
controllers. Based on Figure 7, compared to the selfish controller search
conflicts happen less when the cooperative controller is used. Thus each time
the supervisory controller is triggered, by improving the global performance
of the SaR system, it reduces the number of future search conflicts.

Based on the results and discussions given above, the cooperative controller
significantly outperforms the other methods. The next best controller is the
selfish controller, that is the decentralised control system that remains when
the supervisory MPC controller is excluded. These results further confirm the
robustness of the proposed control architecture to failure of the supervisory
MPC controller, i.e., while the performance degrades after the supervisory MPC
controller is excluded from the control architecture, the performance of the
SaR robotic team is still better than the other control methods used in the
case study.

###  4.3 Structured Simulation Scenarios

Table 6: Victim detection results for the combined scenario (case 5).

|  | Selfish controller |  |  |  |  | Cooperative controller |
---|---|---|---|---|---|---|---|---
victim | health state | number of visits | detection time step |  | victim | health state | number of visits | detection time step
$\boldsymbol{\mathrm{v}_{1}}$ | 6.94 | 0 | - |  | $\boldsymbol{\mathrm{v}_{1}}$ | 6.94 | 0 | -
$\boldsymbol{\mathrm{v}_{2}}$ | 24.9 | 2 | 10 |  | $\boldsymbol{\mathrm{v}_{2}}$ | 17.2 | 1 | 22
$\boldsymbol{V_{3}}$ | 17.97 | 2 | 3 |  | $\boldsymbol{V_{3}}$ | 14.34 | 2 | 8
$\boldsymbol{V_{4}}$ | 0 | 0 | - |  | $\boldsymbol{V_{4}}$ | 6.03 | 1 | 11
$\boldsymbol{V_{5}}$ | 6.94 | 0 | - |  | $\boldsymbol{V_{5}}$ | 21.87 | 1 | 15
$\boldsymbol{V_{6}}$ | 6.94 | 0 | - |  | $\boldsymbol{V_{6}}$ | 11.15 | 1 | 30
|  |  |  |  |  |  |  |

Figure 8: Case 1: Conflict in victim detection.

Figure 9: Case 2: Conflict in victim detection and area coverage.

Figure 10: Case 3: Exploitation of perception fields.

Figure 11: Case 4: Exploitation of sensor accuracies.

Figure 12: Case 5: Combined scenario.

Figure 13: Case 1 - Path taken by SaR robots using the selfish controller
(left) and the cooperative controller (right).

Figure 14: Case 2 - Path taken by SaR robots using the selfish controller
(left) and the cooperative controller (right).

Figure 15: Case 3 - Path taken by SaR robots using the selfish controller
(left) and the cooperative controller (right).

Figure 16: Case 4 - Path taken by SaR robots using the selfish controller
(left) and the cooperative controller (right).

Figure 17: Case 5 - Path taken by SaR robots using the selfish controller
(top) and the cooperative controller (bottom).

Figure 18: Case 1 - Evolution of the scan certainty in time.

Figure 19: Case 2 - Evolution of the scan certainty in time.

Figure 20: Case 3 - Evolution of the scan certainty in time.

Figure 21: Case 4 - Evolution of the scan certainty in time.

Figure 22: Case 5 - Evolution of the scan certainty in time.

In order to further assess the performance of the cooperative controller in a
more structured way and to assess the problem solving behaviour of the
controller when several types of conflicts among the local controllers exist,
five cases of conflicts in smaller scales than the previous simulation
scenarios are considered (see Figures 12-12).

Case 1. Conflict in victim detection: Consider $2$ SaR robots and $2$ victims
in a partially known environment of size $15\times 15$ with $c(x,y,0)=0.2$ for
all $(x,y)\in E$ (see Figure 12). The health states of the victims are $10$
and $50$. Figures 15 and 22 show, respectively, the paths taken by the robots
for $9$ time steps using the cooperative and the selfish control methods, and
the change in the total scan certainty in time.

Based on Figure 15 with the selfish controller, both robots prioritise victim
$\mathrm{v}_{1}$ over victim $\mathrm{v}_{2}$, since these robots are steered
by local controllers that follow (9)-(12), which prioritise visiting a cell
that includes a victim with the worst health state and that is closer to the
SaR robot (where the importance of each factor depends on the values for
parameters $c_{1}$ and $c_{2}$). Since the cell that embeds victim
$\mathrm{v}_{1}$ meets both conditions, victim $\mathrm{v}_{1}$ is the target
of both SaR robots and is detected by them at time step $3$. Afterwards, the
robots continue exploring the environment without moving to victim
$\mathrm{v}_{2}$, who remains outside of their perception fields.

With the cooperative controller, however, SaR robot $\textrm{a}_{2}$ detects
victim $\mathrm{v}_{1}$ at time step $3$ and SaR robot $\textrm{a}_{1}$
detects victim $\mathrm{v}_{2}$ at time step $5$ (see Figures 12 and 15). This
shows that the supervisory MPC controller has successfully coordinated the
actions of the SaR robots in favour of detecting more victims within a given
time span. More specifically, the second term in the objective function of the
MPC controller (see (3.2)) prevents the two SaR robots to cover the same sub-
area of the environment. Since the local loss (considered by the first term in
(3.2)) for redirecting SaR robot $\textrm{a}_{1}$ towards victim
$\textrm{v}_{2}$ is less than that of redirecting the other robot, the
supervisory MPC controller changes the path that has been proposed by the
local controller of SaR robot $\textrm{a}_{1}$ in order to achieve a higher
global gain considering the scan certainty of the environment. The overall
scan certainty for the cooperative controller is $17.2\%$ larger than that of
the selfish controller (see Figure 22).

Case 2. Conflict in victim detection and area coverage: Consider $2$ SaR
robots and $1$ victim in an environment of size $15\times 15$ that is known
except for sub-area $E^{\prime}\subset E$ that is completely unknown, i.e.,
$c(x,y,0)=0$ for all $(x,y)\in E^{\prime}$ and $c(x,y,0)=1$ for all $(x,y)\in
E\setminus E^{\prime}$ (see Figure 12). Figures 15 and 22 show, respectively,
the paths taken by the SaR robots for $9$ time steps using the cooperative and
selfish controllers and the change in the total scan certainty in time.

Figure 15 shows that for the selfish controller, both robots move towards
victim $\mathrm{v}_{1}$ and visit the victim at time step $4$. Afterwards,
both robots follow identical paths to approach the unknown sub-area. These
individual behaviours are steered by the local controllers that, according to
(9)-(12), enforce each robot to prioritise cells that may embed a victim, are
closer to the robot, and gain a higher percentage for the overall scan
certainty of the environment. Since the cell that embeds the victim, in
addition is closer to the robots than the cells within sub-area $E^{\prime}$,
this cell for the local controllers has a higher priority. After visiting the
victim, based on (9), visiting the closest cell within sub-area $E^{\prime}$
becomes the priority of the local controllers.

With the cooperative controller SaR robot $\textrm{a}_{1}$ that is farther
from the victim is redirected via the supervisory MPC controller to move
directly towards sub-area $E^{\prime}$, while SaR robot $\textrm{a}_{2}$ moves
towards victim $\mathrm{v}_{1}$ and visits the victim at time step $4$. With
both the selfish and the cooperative controllers, victim $\mathrm{v}_{1}$ is
detected equally fast, while based on Figure 22, with the cooperative
controller the overall scan certainty is almost $4$ times larger than the
value for the selfish controller. This is an influence of including the second
term of (3.2) in the objective function of the supervisory MPC controller.

Case 3. Exploitation of perception fields: Consider $2$ SaR robots and $2$
victims in a partially known environment of size $15\times 15$ with
$c(x,y,0)=0.5$ for all $(x,y)\in E$ and a set of obstacles shown in black in
Figure 12. SaR robot $\textrm{a}_{1}$ has a sufficiently large perception
field such that it detects both victims, whereas SaR robot $\textrm{a}_{2}$
detects victim $\mathrm{v}_{1}$ only. The health states of victims
$\mathrm{v}_{1}$ and $\mathrm{v}_{2}$ are, respectively, $20$ and $15$. Note
that it is assumed that although the cells that include obstacles (shown in
black) block the movement of the SaR robot, but they do not obstruct the view
of the robot. Figures 15 and 22 show, respectively, the paths taken by the
robots for $10$ time steps using the cooperative and selfish controllers and
the change in the total scan certainty in time.

With the selfish controller, both SaR robots visit victim $\mathrm{v}_{1}$,
which occurs at time steps $2$ and $5$ for robots $\textrm{a}_{1}$ and
$\textrm{a}_{2}$, respectively (see Figure 15). Although the health state of
victim $\mathrm{v}_{1}$ is less critical compared to that of victim
$\mathrm{v}_{2}$, victim $\mathrm{v}_{1}$ is prioritised over victim
$\mathrm{v}_{2}$ by the local controller of SaR robot $\textrm{a}_{2}$,
because this victim can be reached faster and thus the corresponding path
receives a larger grade using (12). Afterwards, the SaR robots continue
exploring the environment without detecting victim $\mathrm{v}_{2}$ due to
their limited perception fields.

With the cooperative controller, SaR robots $\textrm{a}_{1}$ and
$\textrm{a}_{2}$ find victims $\mathrm{v}_{2}$ and $\mathrm{v}_{1}$ at time
steps $10$ and $5$, respectively. Although SaR robot $\textrm{a}_{1}$ is able
to reach victim $\mathrm{v}_{1}$ in a shorter time, the global decision of the
supervisory controller allows victim $\mathrm{v}_{1}$ to be detected later in
order to make sure that victim $\mathrm{v}_{2}$ is detected in time (this is
taken care of via the first term in (3.2), which considers the global gain of
the path grades, instead of the individual/local ones). Thus this simulation
highlights the ability of the cooperative controller to determine locally sub-
optimal tasks for SaR robots, in order to maximise the global mission
performance. Figure 22 shows that by spreading out the SaR robots over the
environment, the cooperative controller achieves an overall scan certainty
that is $18.7\%$ larger than that of the selfish controller.

Case 4: Exploitation of sensor accuracies: Consider $2$ SaR robots in an
environment of size $15\times 15$, where the robots should scan a partially
known environment with two sub-areas $E^{\prime}_{1},E^{\prime}_{2}\subset E$
(see Figure 12). The scan certainty at the initial time step for sub-areas
$E^{\prime}_{1}$ and $E^{\prime}_{2}$ is $0$ and $0.3$, respectively, and for
all cells of $E$ outside these two sub-areas is $0.9$. Figures 17 and 22 show,
respectively, the paths taken by the robots for $10$ time steps using both
cooperative and selfish controllers and the change in the total scan certainty
in time.

Based on Figure 17, with the selfish controller both SaR robots move to sub-
area $E^{\prime}_{1}$ to yield the highest gain in the scan certainty. With
the cooperative controller, however, the robots move to sub-areas
$E^{\prime}_{1}$ and $E^{\prime}_{2}$. Since SaR robot $\textrm{a}_{1}$ is
closer to sub-area $E^{\prime}_{1}$ and has a higher sensor accuracy, it is
sent to sub-area $E^{\prime}_{1}$ by the cooperative controller to yield a
larger overall scan certainty. Based on Figure 22 the overall scan certainty
of the cooperative controller is $29.8\%$ larger than that of the selfish
controller.

Case 5. Combined scenario: Consider $2$ SaR robots, $6$ victims, and a set of
obstacles in an environment of size $30\times 15$, where the robots should
scan a partially known environment with $7$ sub-areas
$E^{\prime}_{1},\ldots,E^{\prime}_{7}\subset E$ (see Figure 12). The scan
certainties at the initial time step for sub-area $E^{\prime}_{1}$ is $0.7$,
for sub-area $E^{\prime}_{2}$ is $0.4$, for sub-areas
$E^{\prime}_{3},E^{\prime}_{4},E^{\prime}_{5},E^{\prime}_{6}$ is $0.2$, for
sub-area $E^{\prime}_{7}$ is $0.1$, and for all cells of $E$ outside these
seven sub-areas is $0.5$. The selfish and cooperative paths of both SaR robots
are shown in Figure 17 for $35$ time steps. Additionally, Table 6 shows the
health state of the victims at the end of the simulation, the number of times
a victim has been visited, and the time step when each victim was first
detected. The change in total scan certainty in time is also illustrated in
Figure 22.

Figure 17 shows that with the selfish controller, both SaR robots visit victim
$\mathrm{v}_{3}$, and then move towards the southwest quadrant of the SaR
environment, where they individually visit victim $\mathrm{v}_{2}$. While $2$
victims are visited doubly by the robots, $4$ victims remain undetected and
$1$ victim deceases. With the cooperative controller, SaR robots
$\textrm{a}_{1}$ and $\textrm{a}_{2}$ visit victims $\mathrm{v}_{4}$ and
$\mathrm{v}_{3}$, respectively. Next they explore different sub-areas of the
environment and detect additional victims $\mathrm{v}_{2}$, $\mathrm{v}_{5}$,
and $\mathrm{v}_{6}$. At the end only $1$ victim remains undetected and no
victim is deceased.

The SaR system detects more victims with the cooperative controller, and no
victim is visited twice, implying the victim search efficiency. Moreover, the
overall scan certainty for the cooperative controller (see Figure 22) is
$27.6\%$ larger than that for the selfish controller.

##  5 Conclusions and Topics for Future Research

Autonomous multi-robot systems are expected to map unknown search-and-rescue
(SaR) environments in a fast and effective way. We have introduced a novel
approach for coordinated mission planning of multi-robot systems for multi-
objective (combined coverage and target-oriented) SaR. The developed control
approach effectively incorporates non-homogeneous imperfect perception
capabilities of the sensors of different robots in order to improve their
performance with respect to the victim detection and area coverage.

The key contributions of the paper are two-fold: in multi-agent control
systems and in search-and-rescue (SaR) robotics. From the point-of-view of
multi-agent control systems, we propose a novel control architecture and
formulation that exploit the imperfect perception capabilities of agents,
coordinate their decisions, and provide a balanced trade-off among various
control objectives in a computationally efficient way. As is also supported by
our simulation results, the developed control system benefits from both
computational efficiency of decentralised control methods and global vision of
centralised control approaches. Additionally, the supervisory level improves
the global control performance based on a predictive and optimal computation
scheme, while local controllers independently steer the agents. Therefore,
although the performance will expectedly degrade, the function of the multi-
agent control system is robust to the failure of this centralised controller.
Furthermore, the integrated formulation proposed in this research allows to
incorporate both expert knowledge (via the fuzzy logic control systems) and
the optimality and predictive capabilities of model predictive control (MPC)
into the decision making of autonomous robots. These contributions are
significant for SaR applications, because existing control methods are mainly
focused on either coverage or target-oriented SaR. Moreover, MPC, which is a
precise control method that systematically handles state and input constraints
and that can provide robustness to SaR uncertainties, has been ignored in the
literature for the crucial task of area coverage in SaR. Our novel approach
and formulation for multi-agent control systems enables MPC to provide all its
strong points for, not only target-oriented, but also coverage-oriented SaR.

We have compared the performance of the resulting cooperative control system
with those of a decentralised selfish control system that excludes the MPC
controller, a pure MPC controller, an ant-colony-based controller, and an
exhaustive random search controller. In $20$ simulated scenarios with randomly
positioned obstacles and victims, the hierarchical control approach showed the
best performance in terms of victim detection efficiency and area coverage.
Moreover, $5$ structured scenarios were designed to simulate conflicting
scenarios and to illustrate the importance of the proposed mathematical
formulations in application. The results proved that in case of conflicts, the
proposed hierarchical controller significantly outperforms the decentralised
controller with a comparable computation time. Moreover, the hierarchical
controller successfully exploits the non-homogeneous perception capabilities
of robots, which improves the overall performance.

In the future, more detailed models that consider the behaviour, physical
capabilities, and intentions of victims for their movement patterns can be
considered. Furthermore, a systematic discussion and evaluation of the
robustness of the proposed control approaches with respect to various sources
of uncertainties, especially uncertainties in the movement of victims, is a
topic of interest for future research. Moreover, in addition to non-
homogeneous perception capabilities, differences in the speed, degrees of
freedom, computational capacity, tasks, and maneuverability of search-and-
rescue robots can be considered. Additionally, combining autonomous learning
methods within the proposed architecture is an interesting topic for future
research. While using a learning-based approach alone may correspond to some
risks for search-and-rescue applications, including such algorithms in a
combined framework, similar to the one proposed in this paper, can result in a
promising performance with adaptability capabilities. Finally, in real-life
implementations the large size of SaR environments increases the computational
burden of the supervisory MPC controller. To address this issue and also to
mitigate the risk of performance degradation due to failure of the supervisory
control level, a similar control architecture with more levels of control may
be proposed. Thus, between the supervisory control level and the steering
local control level, extra levels of control with several distributed MPC
controllers are considered, where each MPC controller supervises a combination
of local sub-areas.

##  6 Statements and Declarations

#### Funding

This research has been supported by the NWO Talent Programme Veni project
“Autonomous drones flocking for search-and-rescue” (18120), which has been
financed by the Netherlands Organisation for Scientific Research (NWO).

#### Competing Interests

The authors declare that they have no financial or non-financial conflict of
interest.

#### Data Availability

The data points, files, and codes for creating the figures represented in the
results of this article are available online at
https://figshare.com/s/9762330a7473363433ab.

#### Author Contributions

Author C. de Koning contributed to designing and implementing the experiments.
Authors C. de Koning and A. Jamshidnejad contributed to the analysis and
interpretation of the results, development of the theoretical contributions,
and composition of the manuscript. Author C. de Koning prepared the first
draft of the manuscript. Author A. Jamshidnejad supervised the study design,
has critically reviewed and edited the manuscript and has prepared the final
version of the paper. Both authors have approved the final version of the
manuscript.

#### Ethics approval

Not applicable.

#### Consent to participate

Not applicable.

#### Consent for publication

Not applicable.

## References

  *   * (1) Casper, J., Murphy, R.R.: Human-robot interactions during the robot-assisted urban search and rescue response at the world trade center. IEEE Transactions on Systems, Man, and Cybernetics, Part B (Cybernetics) 33(3), 367–385 (2003)
  * (2) Coburn, A.W., Spence, R.J.S., Pomonis, A.: Factors determining human casualty levels in earthquakes: mortality prediction in building collapse. In: Proceedings of the Tenth World Conference on Earthquake Engineering, vol. 10, pp. 5989–5994 (1992). Rotterdam, Netherlands
  * (3) Riley, J.M., Endsley, M.R.: The hunt for situation awareness: Human-robot interaction in search and rescue. In: Proceedings of the Human Factors and Ergonomics Society Annual Meeting, vol. 48, pp. 693–697 (2004). Los Angeles, CA
  * (4) Shimanski, C.: Situational awareness in search and rescue operations. In: International Technical Rescue Symposium (2005)
  * (5) Chandarana, M., Hughes, D., Lewis, M., Sycara, K., Scherer, S.: Planning and monitoring multi-job type swarm search and service missions. Journal of Intelligent & Robotic Systems 101(44), 1–14 (2021)
  * (6) Jamshidnejad, A., Frazzoli, E.: Adaptive optimal receding-horizon robot navigation via short-term policy development. In: 15th International Conference on Control, Automation, Robotics and Vision, pp. 21–28. IEEE, Singapore (2018)
  * (7) Beck, Z., Teacy, W.L.T., Jennings, N.R., Rogers, A.C.: Online planning for collaborative search and rescue by heterogeneous robot teams. In: Proceedings of the International Conference on Autonomous Agents & Multiagent Systems, pp. 1024–1033 (2016). Singapore
  * (8) de Alcantara Andrade, F.A., Hovenburg, A.R., de Lima, L.N., Rodin, C.D., Johansen, T.A., Storvold, R., Correia, C.A.M., Haddad, D.B.: Autonomous unmanned aerial vehicles in search and rescue missions using real-time cooperative model predictive control. Sensors 19(19), 4067 (2019)
  * (9) San Juan, V., Santos, M., Andújar, J.M.: Intelligent UAV map generation and discrete path planning for search and rescue operations. Complexity 2018 (2018)
  * (10) Yao, P., Zhao, Z.: Improved glasius bio-inspired neural network for target search by multi-agents. Information Sciences 568, 40–53 (2021)
  * (11) Galceran, E., Carreras, M.: A survey on coverage path planning for robotics. Robotics and Autonomous systems 61(12), 1258–1276 (2013)
  * (12) Koenig, S., Liu, Y.: Terrain coverage with ant robots: a simulation study. In: Proceedings of the Fifth International Conference on Autonomous Agents, pp. 600–607 (2001). Montreal, Canada
  * (13) Wagner, I.A., Altshuler, Y., Yanovski, V., Bruckstein, A.M.: Cooperative cleaners: A study in ant robotics. The International Journal of Robotics Research 27(1), 127–151 (2008)
  * (14) Yang, S.X., Luo, C.: A neural network approach to complete coverage path planning. IEEE Transactions on Systems, Man, and Cybernetics, Part B (Cybernetics) 34(1), 718–724 (2004)
  * (15) Yang, Y., Polycarpou, M.M., Minai, A.A.: Multi-UAV cooperative search using an opportunistic learning method. Journal of Dynamic Systems, Measurement, and Control 129 (2007)
  * (16) Tutsoy, O., Barkana, D.E., Balikci, K.: A novel exploration-exploitation-based adaptive law for intelligent model-free control approaches. IEEE Transactions on Cybernetics, 1–9 (2021). https://doi.org/10.1109/TCYB.2021.3091680
  * (17) Tutsoy, O., Barkana, D.E., Colak, S.: Learning to balance a Nao robot using reinforcement learning with symbolic inverse kinematic. Transactions of the Institute of Measurement and Control 39(11), 1735–1748 (2017)
  * (18) Arnold, R., Jablonski, J., Abruzzo, B., Mezzacappa, E.: Heterogeneous UAV multi-role swarming behaviors for search and rescue. In: IEEE Conference on Cognitive and Computational Aspects of Situation Management, pp. 122–128 (2020). Victoria, BC, Canada
  * (19) Farrokhsiar, M., Pavlik, G., Najjaran, H.: An integrated robust probing motion planning and control scheme: A tube-based MPC approach. Robotics and Autonomous Systems 61(12), 1379–1391 (2013)
  * (20) Hoy, M., Matveev, A.S., Savkin, A.V.: Collision free cooperative navigation of multiple wheeled robots in unknown cluttered environments. Robotics and Autonomous Systems 60(10), 1253–1266 (2012)
  * (21) Cooper, J.R.: Optimal multi-agent search and rescue using potential field theory. In: AIAA Scitech 2020 Forum, p. 0879 (2020)
  * (22) Paez, D., Romero, J.P., Noriega, B., Cardona, G.A., Calderon, J.M.: Distributed particle swarm optimization for multi-robot system in search and rescue operations. IFAC-PapersOnLine 54(4), 1–6 (2021)
  * (23) Choi, H., Brunet, L., How, J.P.: Consensus-based decentralized auctions for robust task allocation. IEEE Transactions on Robotics 25(4), 912–926 (2009)
  * (24) Liu, Y., Nejat, G.: Robotic urban search and rescue: A survey from the control perspective. Journal of Intelligent & Robotic Systems 72(2), 147–165 (2013)
  * (25) Best, G., Hollinger, G.A.: Decentralised self-organising maps for multi-robot information gathering. In: 2020 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS), pp. 4790–4797. IEEE/RSJ, Las Vegas (2020)
  * (26) Otte, M., Kuhlman, M.J., Sofge, D.: Auctions for multi-robot task allocation in communication limited environments. Autonomous Robots 44, 547–584 (2020)
  * (27) Kashino, Z., Nejat, G., Benhabib, B.: Aerial wilderness search and rescue with ground support. Journal of Intelligent & Robotic Systems 99, 147–163 (2020)
  * (28) Tol, D., Hoekstra, J., Jamshidnejad, A.: A bi-level local and global model predictive control architecture for air traffic management. In: 2021 IEEE International Intelligent Transportation Systems Conference (ITSC), pp. 361–365 (2021)
  * (29) Khamis, A.M., Elmogy, A.M., Karray, F.O.: Complex task allocation in mobile surveillance systems. Journal of Intelligent & Robotic Systems 64(1), 33–55 (2011)
  * (30) Elston, J., Frew, E.W.: Hierarchical distributed control for search and tracking by heterogeneous aerial robot networks. In: IEEE International Conference on Robotics and Automation, pp. 170–175 (2008). Pasadena, CA, USA
  * (31) Chandler, P.R., Pachter, M., Rasmussen, S.: UAV cooperative control. In: Proceedings of the American Control Conference, vol. 1, pp. 50–55 (2001). Arlington, VA, USA
  * (32) Grogan, S., Pellerin, R., Gamache, M.: The use of unmanned aerial vehicles and drones in search and rescue operations - a survey. In: Proceedings of the PROLOG (2018). Hull, UK
  * (33) Krzysiak, R., Butail, S.: Information-based control of robots in search-and-rescue missions with human prior knowledge. IEEE Transactions on Human-Machine Systems 52(1), 52–63 (2021)
  * (34) Ganesan, S., Shakya, M., Aqueel, A.F., Nambiar, L.M.: Small disaster relief robots with swarm intelligence routing. In: Proceedings of the 1st International Conference on Wireless Technologies for Humanitarian Relief, pp. 123–127 (2011). Kollam, India
  * (35) Wang, W., Joshi, R., Kulkarni, A., Leong, W.K., Leong, B.: Feasibility study of mobile phone WiFi detection in aerial search and rescue operations. In: Proceedings of the 4th Asia-Pacific Workshop on Systems, pp. 1–6 (2013). Singapore
  * (36) Dousai, N.M.K., Lončarić, S.: Detecting humans in search and rescue operations based on ensemble learning. IEEE Access 10, 26481–26492 (2022)
  * (37) Llasag, R., Marcillo, D., Grilo, C., Silva, C.: Human detection for search and rescue applications with UAVs and mixed reality interfaces. In: 2019 14th Iberian Conference on Information Systems and Technologies (CISTI), pp. 1–6 (2019)
  * (38) Pinheiro, G.P.M., Miranda, R.K., Praciano, B.J.G., Santos, G.A., Mendonça, F.L.L., Javidi, E., da Costa, J.P.J., de Sousa, R.T.J.: Multi-sensor wearable health device framework for real-time monitoring of elderly patients using a mobile application and high-resolution parameter estimation. Frontiers in Human Neuroscience (2022)
  * (39) Hart, P.E., Nilsson, N.J., Raphael, B.: A formal basis for the heuristic determination of minimum cost paths. IEEE transactions on Systems Science and Cybernetics 4(2), 100–107 (1968)
  * (40) Yen, J.Y.: Finding the k shortest loopless paths in a network. Management Science 17(11), 712–716 (1971)
  * (41) Jamshidnejad, A., Papamichail, I., Papageorgiou, M., De Schutter, B.: Sustainable model-predictive control in urban traffic networks: Efficient solution based on general smoothening methods. IEEE Transactions on Control Systems Technology 26(3), 813–827 (2018)
  * (42) Diehl, M., Bock, H.G., Schlöder, J.P.: A real-time iteration scheme for nonlinear optimization in optimal feedback control. SIAM Journal on Control and Optimization 43(5), 1714–1736 (2005)
  * (43) Kreinovich, V., Kosheleva, O., Shahbazova, S.N.: Why triangular and trapezoid membership functions: A simple explanation. In: Shahbazova, S.N., Sugeno, M., Kacprzyk, J. (eds.) Recent Developments in Fuzzy Logic and Fuzzy Sets: Dedicated to Lotfi A. Zadeh, pp. 25–31. Springer, Cham (2020)


# DL-Corrector-Remapper: A grid-free bias-correction deep learning methodology
for data-driven high-resolution global weather forecasting

Tao Ge
Washington University in St. Louis
St. Louis, MO. 63130
<EMAIL_ADDRESS>
&Jaideep Pathak
NVIDIA
Santa Clara, CA. 95050
<EMAIL_ADDRESS>
Akshay Subramaniam
NVIDIA
Santa Clara, CA. 95050
<EMAIL_ADDRESS>
&Karthik Kashinath
NVIDIA
Santa Clara, CA. 95050
<EMAIL_ADDRESS>
Worked performed during an internship at NVIDIA

###### Abstract

Data-driven models, such as FourCastNet (FCN), have shown exemplary
performance in high-resolution global weather forecasting. This performance,
however, is based on supervision on mesh-gridded weather data without the
utilization of raw climate observational data, the gold standard ground truth.
In this work we develop a methodology to correct, remap, and fine-tune gridded
uniform forecasts of FCN so it can be directly compared against observational
ground truth, which is sparse and non-uniform in space and time. This is akin
to bias-correction and post-processing of numerical weather prediction (NWP),
a routine operation at meteorological and weather forecasting centers across
the globe. The Adaptive Fourier Neural Operator (AFNO) architecture is used as
the backbone to learn continuous representations of the atmosphere. The
spatially and temporally non-uniform output is evaluated by the non-uniform
discrete inverse Fourier transform (NUIDFT) given the output query locations.
We call this network the Deep-Learning-Corrector-Remapper (DLCR). The
improvement in DLCR’s performance against the gold standard ground truth over
the baseline’s performance shows its potential to correct, remap, and fine-
tune the mesh-gridded forecasts under the supervision of observations.

##  1 Introduction

Reliable weather forecasting models play a crucial role in preparing for the
harsh consequences of climate change by providing early warnings of extreme
weather events for disaster mitigation. Numerous deep learning (DL)-based
weather prediction models have been developed to forecast global weather under
the supervision of the reanalysis mesh-gridded data [1, 2, 3, 4, 5]. FCN [6,
7] has shown exemplary performance in high-resolution data-driven global
weather forecasting measured by the similarity of its forecasts against mesh-
gridded reanalysis data, ERA5[8]. However, FCN’s mesh-gridded forecasts could
not be directly compared against the raw, sparse, non-uniform observations.
Yet, the gold standard for weather and climate model performance is their
ability to match observations. Further, because FCN is trained on reanalysis
data, it is likely to have biases with respect to direct observations. To
address these issues, we develop a DL-based _grid-free_ model that remaps and
corrects the mesh-gridded forecasts to arbitrary locations in space and time.
This methodology is sufficiently general that it can be applied to _any_
gridded forecast, not just from FCN.

Several data-driven approaches have been developed to refine mesh-gridded
forecasts using sparse, non-uniform observational data [9, 10, 11]. However,
our objective is different from existing work in the following distinct ways:
(1) The desired model should accept arbitrary output query locations, _i.e._
it should be grid-free (existing works are only suitable for fixed locations);
(2) Some observations may be missing randomly over space and time, so the
desired model should handle random missing data seamlessly to fully utilize
the available data at all locations and time steps; (3) The desired model
should be capable of processing large spatio-temporal data from long
forecasts, including all relevant variables, at high resolution, which can be
highly memory-intensive.

To address these challenges, we develop a Fourier-based approach to learn
continuous representations of the input data under the supervision of sparse,
non-uniform observations. The model is grid-free since the continuous data can
be reconstructed from the Fourier coefficients at any location. Moreover,
because the network itself is independent of the output query, the model can
naturally process missing data. To the best of our knowledge, this is the
first time that a _grid-free_ DL model has been developed to re-map and
correct mesh-gridded high-resolution weather forecasts to sparse, non-uniform
observations.

##  2 Methodology and Data

###  2.1 Type II Non-Uniform DFT

We use Fourier coefficients for continuous representations of spatio-temporal
data. Therefore, the projection of the continuous representation to spatial
positions is the real part of type II non-uniform inverse discrete Fourier
transform (NUIDFT), which can be written in matrix form as:

$\displaystyle\cos(2\pi QM^{T})F_{real}-\sin(2\pi QM^{T})F_{img},$ (1)

where $Q\in\mathbb{R}^{N\\!\times\\!2}$ denotes the query matrix,
$M\in\mathbb{R}^{W\\!H\\!\times\\!2}$ is the frequency basis, $F_{real}$ and
$F_{img}\in\mathbb{R}^{W\\!H\\!\times\\!1}$ denote the real and imaginary
parts of the Fourier coefficients, respectively.

The sparse NUIDFT is developed to reduce the memory footprint by dividing the
queries into different groups with respect to their actual length and
truncating zero elements. For very high-resolution data, query locations can
be randomly sub-sampled in each epoch to further reduce memory usage.

###  2.2 DL Corrector-Remapper and Baselines

Figure 1: Flowchart of DL Corrector-Remapper (DLCR). DLCR consists of three
parts: (i) image patching and coding; (ii) U-shaped AFNO blocks; and (iii)
decoding and NUIDFT. The typical size of the output tensor is denoted under
the name of each block.

We develop DL-Corrector-Remapper (DLCR), a Fourier-based network that
incorporates NUIDFT into the backbone of Adaptive Fourier Neural Operators
(AFNO) [12], to remap and bias-correct mesh-gridded forecasts to arbitrary
locations in space and time, under the supervision of the gold standard –
sparse, non-uniform observational data.

The overall structure of the proposed network is shown in Figure 1.
Essentially, the output of the model is the Fourier coefficient, so one can
use arbitrary queries to get arbitrary spatial outputs. Applying the inverse
fast Fourier transform (IFFT) yields mesh-gridded output, while applying the
NUIDFT with the query yields sparse observational output (see bottom right in
Figure 1).

In the _patching process_ , non-overlapping patches of the image are stacked
in the channel dimension to form a token sequence [12]. Then, the channel
coder compresses the channel size to reduce the memory footprint by accounting
for the spatio-temporal and inter-variable dependencies. The _U-shaped AFNO
blocks_ are the backbone of the proposed model. The channel dimension is
halved for each subsequent block $k$ if $k<(K+1)//2$, and doubled if
$k\geq(K+1)//2$, where $K$ denotes the total number of blocks. This U-shaped
structure reduces the memory footprint without sacrificing the depth of the
network. Moreover, _skip connections_ allow feature maps with more channels to
flow in the network, so the narrowest block in the network does not reduce the
complexity of the output.

The proposed network (DLCR) is compared to the NUIDFT-based interpolation and
the widely used standard U-Net. Figure 2 shows a schematic of these models and
the corresponding number of trainable parameters.

Figure 2: Model comparison of the interpolation, U-Net, and DLCR. These three
methods utilize NUIDFT to reconstruct the spatio-temporal data from the
continuous representation. The interpolation method has no trainable
parameters, while U-Net and FCR both have $\sim$3 billion trainable
parameters.

###  2.3 Loss function

To correct biases of mesh-gridded forecasts whilst preserving their spatio-
temporal structure, shapes, and patterns, we develop a novel loss function
with two terms: (i) a magnitude difference term; and (ii) a structural
similarity term,

$\displaystyle\operatorname*{argmin}_{\theta}||$ $\displaystyle
y_{sp}-\Gamma_{\theta}(x,q_{sp})||^{2}_{2}+\lambda(1-LCC(y_{gd},\Gamma_{\theta}(x,q_{gd})))$
(2)

where $y_{sp}$ denotes sparse, non-uniform observational data, $q_{sp}$ is the
corresponding sparse output query, $y_{gd}$ denotes the reanalysis mesh-
gridded data (ERA5), $q_{gd}$ is the mesh-gridded query, $\Gamma_{\theta}$
denotes the trainable model, $x$ denotes the input data, $\lambda$ denotes a
scalar that controls the weight of the structural similarity term, and $LCC$
denotes the local cross-correlation [13] (see Appendix A).

###  2.4 Datasets and Training

The input dataset is FCN inference data from 2000 to 2018 at 0.25°resolution
[6]. We use the global observational data acquired from National Centers for
Environmental Information as the ground truth reference [14]. Four overlapping
variables are selected: surface wind velocity components (U and V), surface
temperature (t2m), and mean sea level pressure (MSLP). A data sample is shown
in Figure 3. The input and reference data are normalized by the global mean
and variance of each variable.

Figure 3: The mesh-gridded input (top) and non-uniform reference (bottom)
data. Wind velocity V is omitted in this figure. The reference data is
extremely sparse and missing randomly.

Each forecast is five-days-long, with a 6-hourly time step, generated across
19 years, which yields 27360 samples in total. We use 23040 samples from 2000
to 2015 to train the model and 1440 samples from the year 2017 for validation.
Moreover, 2% of the observational sites are removed from the training set, so
these positions are ‘unobserved’ or _hidden_. ‘Unobserved’ time is used to
test the reliability of future forecasts, while ‘unobserved’ positions are
used to test the model’s ability to produce observation-quality forecasts for
locations that it has never seen before. Since the reference data is extremely
sparse, the output query is perturbed by Gaussian noise $\mathcal{N}(0,0.04)$
to address spatial over-fitting.

The networks are trained on 64 A100 NVIDIA GPUs with 80 GB memory on the
Selene supercomputer [15]. Each training epoch takes 3000 seconds. The
optimizer is ADAM, and the start learning rate is $10^{-3}$. We train the
model for 10 epochs with a scheduler that reduces the learning rate by a
factor of 0.7 every 3 epochs.

##  3 Results

Figure 4: MSE of the output against the observation vs. prediction time step
(lower is better). Our proposed model (DLCR) outperforms baselines across
observed and unobserved locations.

Figure 4 shows the plots of mean square error (MSE) that are averaged over 80
instances across the year 2017 (validation). The horizontal axis denotes the
time step, where each step corresponds to 12 hours. For observed positions,
the performance of DLCR is close to the performance of U-Net, and they both
outperform the interpolation baseline. For unobserved positions, DLCR
outperforms the interpolation baseline and U-Net, and it performs better on
more complicated variables (U and V), whereas the performance of the U-Net is
even worse than the performance of the interpolation in estimating t2m and
MSLP. Besides, the difference between the solid and dashed blue curves may
indicate a small systematic bias in the unobserved locations. (see appendix
for additional information)

##  4 Conclusion

We propose a grid-free DL model (DLCR) to correct, remap, and fine-tune
gridded uniform high-resolution global weather forecasts under the supervision
of the gold standard ground truth – observations – which are sparse and non-
uniform in space and time, can be noisy, and have gaps and missing data. DLCR
outperforms baselines even when the output positions are unobserved. In this
paper DLCR is developed for correcting and remapping FourCastNet’s output, but
the method itself can be applied to any forecast data. This work has
implications for improving a broad range of data-driven and traditional NWP
forecasting approaches by making their outputs comparable to and trainable on
observational data. Importantly, DLCR allows arbitrary query locations and is
trained to produce _observation-quality_ forecasts in ‘unobserved’ locations,
therefore it provides a way to improve the reliability of data-driven and NWP
forecasts in locations on the globe and for time instances where observations
do not exist.

## References

  * [1] P. D. Dueben and P. Bauer.  Challenges and design choices for global weather and climate models based on machine learning.  Geoscientific Model Development, 11(10):3999–4009, 2018.
  * [2] Jonathan A. Weyn, Dale R. Durran, and Rich Caruana.  Can Machines Learn to Predict Weather? Using Deep Learning to Predict Gridded 500-hPa Geopotential Height From Historical Weather Data.  Journal of Advances in Modeling Earth Systems, 11(8):2680–2693, 2019.
  * [3] Jonathan A. Weyn, Dale R. Durran, and Rich Caruana.  Improving Data-Driven Global Weather Prediction Using Deep Convolutional Neural Networks on a Cubed Sphere.  Journal of Advances in Modeling Earth Systems, 12(9):e2020MS002109, 2020.  e2020MS002109 10.1029/2020MS002109.
  * [4] Stephan Rasp, Peter D. Dueben, Sebastian Scher, Jonathan A. Weyn, Soukayna Mouatadid, and Nils Thuerey.  WeatherBench: A Benchmark Data Set for Data-Driven Weather Forecasting.  Journal of Advances in Modeling Earth Systems, 12(11):e2020MS002203, 2020.
  * [5] Ryan Keisler.  Forecasting Global Weather with Graph Neural Networks, 2022.
  * [6] Jaideep Pathak, Shashank Subramanian, Peter Harrington, Sanjeev Raja, Ashesh Chattopadhyay, Morteza Mardani, Thorsten Kurth, David Hall, Zongyi Li, Kamyar Azizzadenesheli, Pedram Hassanzadeh, Karthik Kashinath, and Animashree Anandkumar.  FourCastNet: A Global Data-driven High-resolution Weather Model using Adaptive Fourier Neural Operators, 2022.
  * [7] Thorsten Kurth, Shashank Subramanian, Peter Harrington, Jaideep Pathak, Morteza Mardani, David Hall, Andrea Miele, Karthik Kashinath, and Animashree Anandkumar.  FourCastNet: Accelerating Global High-Resolution Weather Forecasting using Adaptive Fourier Neural Operators, 2022.
  * [8] H. Hersbach, Berrisford Bell, B., Biavati P., Horányi G., J. A., Muñoz Sabater, J. Nicolas, Peubey, R. C., Radu, I. Rozum, D. Schepers, A. Simmons, C. Soci, D. Dee, and J-N Thépaut.  ERA5 hourly data on single levels from 1959 to present.  N Copernicus Climate Change Service (C3S) Climate Data Store (CDS).  Accessed: 2022-07.
  * [9] Stephan Rasp and Sebastian Lerch.  Neural networks for postprocessing ensemble weather forecasts.  Monthly Weather Review, 146(11):3885 – 3900, 2018.
  * [10] Federico Amato, Fabian Guignard, Sylvain Robert, and Mikhail Kanevski.  A novel framework for spatio-temporal prediction of environmental data using deep learning.  Scientific Reports, 10(22243), 2020.
  * [11] Hui Zhang, Yaqiang Wang, Dandan Chen, Dian Feng, Xiaoxiong You, and Weichen Wu.  Temperature Forecasting Correction Based on Operational GRAPES-3km Model Using Machine Learning Methods.  Atmosphere, 13(2), 2022.
  * [12] John Guibas, Morteza Mardani, Zongyi Li, Andrew Tao, Anima Anandkumar, and Bryan Catanzaro.  Adaptive Fourier Neural Operators: Efficient Token Mixers for Transformers.  CoRR, abs/2111.13587, 2021.
  * [13] B. B. Avants, C. L. Epstein, M. Grossman, and J. C. Gee.  Symmetric diffeomorphic image registration with cross-correlation: evaluating automated labeling of elderly and neurodegenerative brain.  Med Image Anal., 12(1):26–41, 2008.
  * [14] NOAA National Centers for Environmental Information (2001): Global Surface Hourly.  NOAA National Centers for Environmental Information: https://www.ncei.noaa.gov/metadata/geoportal/rest/metadata/item/gov.noaa.ncdc:C00532/html.  Accessed: 2022-07.
  * [15] NVIDIA Selene: Leadership-Class Supercomputing Infrastructure.  https://www.nvidia.com/en-us/on-demand/session/supercomputing2020-sc2019/.

##  Appendix A Appendix: Local Cross-Correlation

Local cross-correlation (LCC) is a similarity metric that focuses on the
similarity of patterns, structures, and high contrast boundaries between
signals, written as

$\displaystyle
LCC(A,B)=\frac{\Big{(}\sum_{i}\sum_{j}\big{(}A[i,j]-\bar{A}[i,j]\big{)}\big{(}B[i,j]-\bar{B}[i,j]\big{)}\Big{)}^{2}}{\sum_{i}\sum_{j}\Big{(}A[i,j]-\bar{A}[i,j]\Big{)}^{2}\sum_{i}\sum_{j}\Big{(}B[i,j]-\bar{B}[i,j]\Big{)}^{2}},$
(3)

where $A$ and $B$ denote two arbitrary 2-dimensional images, and $\bar{A}$
denotes the local mean that is calculated by the convolution of the image $A$
and an averaging kernel $V$, _i.e._ ,

$\displaystyle\bar{A}[i,j]=\sum_{n}\sum_{m}A[i-n,j-m]V[n,m].$ (4)

In this work, V is chosen to be a $5\times 5$ uniform matrix whose sum is
normalized to 1. The size of the kernel controls the scale of the
detectability. LCC with a large kernel only captures large structures and
ignores subtle changes, while LCC with a small kernel captures subtle
variations yet is more vulnerable to noise.

##  Appendix B Appendix: Additional Figures

Figure 5 shows the scatter plot and the interpolated mesh-gridded image of the
reference observation data, as well as the scatter plot and the interpolated
output inference with a random query (6000 locations). The display region
covers -40°to 40°of latitude and -22.5°to 50°of longitude, which corresponds
to the continent of Africa. These query positions are not necessary to be on
the mesh grid or observed. The randomly queried output of the proposed method
still looks realistic in the region with a limited number of observational
sites.

Figure 5: Figures of the reference and DLCR inference output for the region of
Africa. The DLCR inference output is generated via a random query unseen by
the model. It is clear that interpolation of DLCR’s output, after querying
DLCR densely over the African continent, has significantly better qualitative
performance in regions of Africa with very sparse observations, whereas
interpolation of the raw observations produces unphysical results in these
regions. Thus the benefit of a grid-free model that can be queried at
arbitrary locations is evident.

Figure 6 shows the observation positions in this work, where hidden positions
are labelled as orange triangles. The world map is evenly divided into 48
patches, and these hidden positions are randomly selected in each patch.

Figure 6: Figure showing the locations of the surface stations. Observed
positions are indicated by grey dots, and unobserved positions are indicated
by orange pins.