An opportunistic linear-convex algorithm for localization in mobile robot networks

TL;DR

This paper introduces a distributed, opportunistic linear algorithm for localizing mobile robots in a network, which converges asymptotically and is robust to noise, especially with known beacon nodes.

Contribution

It presents a novel linear-convex, distributed localization algorithm for mobile robot networks that operates opportunistically and guarantees convergence under certain conditions.

Findings

Algorithm converges asymptotically to true robot locations.

Performance is robust to measurement noise with proposed modifications.

Requires at least one beacon for precise tracking of mobile networks.

Abstract

In this paper, we develop a \textcolor{black}{\emph{distributed}} algorithm to localize a network of robots moving arbitrarily in a bounded region. In the case of such mobile networks, the main challenge is that the robots may not be able to find nearby robots to implement a distributed algorithm. We address this issue by providing an opportunistic algorithm that only implements a location update when there are nearby robots and does not update otherwise. We assume that each robot measures a noisy version of its motion and the distances to the nearby robots. To localize a network of mobile robots in~$\mathbb{R}^m$, we provide a simple \emph{linear} update, which is based on barycentric coordinates and is linear-convex. We abstract the corresponding localization algorithm as a Linear Time-Varying (LTV) system and show that it asymptotically converges to the true locations~of~the robots. We first focus on the noiseless case, where the distance and motion vectors are known (measured) perfectly, and provide sufficient conditions on the convergence of the algorithm. We then evaluate the performance of the algorithm in the presence of noise and provide modifications to counter the undesirable effects of noise. \textcolor{black}{We further show that our algorithm precisely tracks a mobile network as long as there is at least one known beacon (a node whose location is perfectly known).

Equations

$${\mathbf{x}}_{k+1}^{i\ast }={\mathbf{x}}_k^{i\ast }+{\tilde{\mathbf{x}}}_{k+1}^i,i\in \Psi ,$$

$${\hat{\mathbf{x}}}_k^i={\tilde{\mathbf{x}}}_k^i+n_k^i,$$

$$\widehat{d}^{ij}_{k}={{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\widetilde{d}^{ij}_{k}}}+r^{ij}_{k},$$

$$i\in \mathcal{C}({\Theta }_i(k)),\text{mbox}if\sum _{j\in {\Theta }_i(k)}{A}_{{\Theta }_i(k)\cup \{i\}\setminus j}=A_{{\Theta }_i(k)},$$

$$i\notin \mathcal{C}({\Theta }_i(k)),\text{mbox}if\sum _{j\in {\Theta }_i(k)}{A}_{{\Theta }_i(k)\cup \{i\}\setminus j}>A_{{\Theta }_i(k)},$$

$${\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}_{k}^{i\ast}=a_{k}^{{i}1}\mathbf{x}_{k}^{1\ast}+a_{k}^{{i}2}\mathbf{x}_{k}^{2\ast}+a_{k}^{{i}3}\mathbf{x}_{k}^{3\ast}=\sum_{j\in\Theta_{i}(k)}a_{k}^{ij}\mathbf{x}_{k}^{j\ast},}$$

$${\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}a_{k}^{{i}j}=\frac{A_{\Theta_{i}(k)\cup\{i\}\setminus j}}{A_{\Theta_{i}(k)}}}.$$

$$\displaystyle\mathbf{x}^{i}_{k+1}=\alpha_{k}\mathbf{x}_{k}^{i}+(1-\alpha_{k})\sum_{j\in\Theta_{i}(k)}{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}a_{k}^{ij}}\mathbf{x}_{k}^{j}+\widetilde{\mathbf{x}}_{k+1}^{i},$$

$$\displaystyle\alpha_{k}=\left\{\begin{array}[]{ll}1,&\forall k~{}|~{}\Theta_{i}(k)=\emptyset,\\ \in\left[\beta,1\right),&\forall k~{}|~{}\Theta_{i}(k)\neq\emptyset,\end{array}\right.$$

$${\mathbf{x}}_{k+1}^i$$

$$\displaystyle{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}a_{k}^{ij}}=\begin{cases}{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}p_{k}^{ij}},&\mbox{if}~{}j\in\Theta_{i}(k)\cap\Omega,\\ {\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}b_{k}^{im}},&\mbox{if}~{}m\in\Theta_{i}(k)\cap\kappa.\end{cases}$$

$$\displaystyle 0<\beta\leq{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}p_{k}^{ii}}\leq 1,\qquad\forall k,i\in\Omega.$$

$${\mathbf{x}}_{k+1}={\mathbf{P}}_k{\mathbf{x}}_k+{\mathbf{B}}_k{\mathbf{u}}_k+{\tilde{\mathbf{x}}}_{k+1},k>0,$$

$$\displaystyle 0<\alpha\leq{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}b_{k}^{im}},\qquad\forall k,i\in\Omega,{m\in\Theta_{i}(k)\cap\kappa},$$

$$\displaystyle{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}b_{k}^{im}}=\frac{A_{\Theta_{i}(k)\cup\{i\}\setminus m}}{A_{\Theta_{i}(k)}},\qquad m\in\kappa,$$

$${\mathbf{x}}_{k+1}={\mathbf{P}}_k{\mathbf{x}}_k,$$

$$0<{\beta }_1\le p_k^{ii},{\beta }_1\in \mathbb{R}.$$

$$\sum _j{p}_k^{ij}\le {\beta }_2<1,{\beta }_2\in \mathbb{R}.$$

$$\prod _t{\mathbf{M}}_t=\prod _k{\mathbf{P}}_k.$$

$${\Vert M_j\Vert }_{\infty }\le 1-{{\beta }_1}^{|M_j|-1}{\beta }_2.$$

$$\displaystyle{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}|M_{j}|\leq L<{\infty},\qquad\forall j,~{}L\in\mathbb{N};}$$

$$|M_j|\le L_1<\infty ,\forall M_j\in J_1,$$

$$|M_j|<\infty ,\forall M_j\notin J_1;$$

$$\exists M_j\in J_2:|M_j|\le \frac{1}{\ln \left({\beta }_1\right)}\ln \left(\frac{1-e^{(-{\gamma }_2{i}^{-{\gamma }_1})}}{1-{\beta }_2}\right)+1,$$

$${\Vert \mathbf{y}(t+1)\Vert }_{\infty }$$

$$\displaystyle\sum_{j\in\Theta_{i}(k)\cap\Omega}{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}p_{k}^{ij}}=1-\sum_{m\in\Theta_{i}(k)\cap\kappa}{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}b_{k}^{im}},$$

$$\displaystyle\sum_{j\in\Theta_{i}(k)\cap\Omega}{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}p_{k}^{ij}}\leq 1-\alpha<1,$$

$${\mathbf{x}}_{k+1}={\mathbf{P}}_k{\mathbf{x}}_k+{\mathbf{B}}_k{\mathbf{u}}_k+{\tilde{\mathbf{x}}}_{k+1},k>0,$$

$$\underset{k\to \infty }{\lim }\prod _{l=0}^k{\mathbf{P}}_l=0_{N\times N}.$$

$${\mathbf{x}}_{k+1}^{\ast }={\mathbf{P}}_k{\mathbf{x}}_k^{\ast }+{\mathbf{B}}_k{\mathbf{u}}_k+{\tilde{\mathbf{x}}}_{k+1}.$$

$${\mathbf{e}}_{k+1}\triangleq {\mathbf{x}}_{k+1}^{\ast }-{\mathbf{x}}_{k+1}={\mathbf{P}}_k({\mathbf{x}}_k^{\ast }-{\mathbf{x}}_k)={\mathbf{P}}_k{\mathbf{e}}_k,$$

Full text

We first focus on the noiseless case, where the distance and motion vectors are known (measured) perfectly, and provide sufficient conditions on the convergence of the algorithm. We then evaluate the performance of the algorithm in the presence of noise and provide modifications to counter the undesirable effects of noise. We further show that our algorithm precisely tracks a mobile network as long as there is at least one known beacon (a node whose location is perfectly known).

I Introduction

In order to navigate reliably and perform useful tasks in robotic networks, a mobile robot must know its exact location. Robot localization, estimation of a robot’s location from sensor data, is thus a fundamental problem in providing autonomous capabilities to a mobile robot. Although some robotic systems rely on a Global Positioning System (GPS) to determine their location in a global reference frame, it is impractical to use GPS in many indoor applications. Localization problems in robotics can be broadly grouped into two categories: position tracking; and global localization. Position tracking, e.g., [1], also known as local localization, requires the knowledge of robots’ starting locations, whereas in global localization, e.g., [2], no prior estimate of the initial locations exists. In terms of the information used for estimating the locations, localization schemes can be classified as range-based,[3, 4], and range-free,[5, 6]. While the former depends on measuring the distance and/or angle between the nodes, the latter makes no assumptions about the availability of such information and relies on the connectivity of the network.

Since the early work on navigation with autonomous mobile robots, a variety of centralized and distributed techniques have been proposed to tackle the localization problem. Despite higher accuracy in small-sized networks, centralized localization schemes suffer from scalability issues, and are not feasible in large-scale networks. In addition, comparing to the distributed algorithms, the centralized schemes are less reliable and require higher computational complexity due to, e.g., the accumulated inaccuracies caused by multi-hop communications over a wireless network.

Robot localization approaches include but are not limited to dead-reckoning, [7, 8, 9], Simultaneous Localization and Mapping (SLAM),[10, 11, 12, 13, 14], Monte Carlo techniques [5, 6, 3, 15, 16], and Kalman Filtering methods [17, 18, 19, 20, 21, 22, 23, 24, 25]; other related works include[26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38]. We briefly describe the related work below.

Dead-reckoning, [7, 8], is a common method to estimate the location of a mobile robot. It uses the wheel rotation measurements to compute the offset from a known starting position. Despite the low cost, simplicity, and easy implementation in real time, dead reckoning methods are prone to accuracy problems due to accumulating wheel slippage errors, which grow without bound over time. Therefore, these methods are only suitable for applications where the robots have good estimates of their initial locations, and their tasks involve exploring only short distances.

When the map of the environment is not available a priori, Simultaneous Localization and Mapping techniques,[10, 11, 12, 13], can be used to build a map of an unexplored environment by a mobile robot, while simultaneously navigating the environment using the map. The main disadvantage of most SLAM algorithms is the high computational complexity, which makes them less efficient specially in larger multi-robot networks. Another disadvantage of SLAM-based solutions is the sharing of the map when multiple robots are involved.

When ranging data is noisy, estimation-based localization techniques are widely used. Sequential Bayesian Estimation (SBE) methods use the recursive Bayes rule to estimate the likelihood of a robot’s location. The solution to SBE is generally intractable and cannot be determined analytically. An alternative approach is Kalman-based techniques, which are only optimal when the uncertainties are Gaussian and the system dynamics are linear. However, localization has always been considered as a nonlinear problem and hence the optimality of Kalman-based solutions are not guaranteed.

To address the nonlinear nature of localization problems, other suboptimal solutions to approximate the optimal Bayesian estimation include Particle Filters (PF) and Extended Kalman Filters (EKF). In particular, Sequential Monte Carlo (SMC) method is a PF that exploits posterior probability to determine the future location of a robot. Monte Carlo Localization (MCL) algorithms can solve global localization in a robust and efficient way. For example, Ref. [5] introduces the Monte Carlo Localization (MCL) method, which exploits mobility to improve the accuracy of localization. Inspired by [5], the authors in [6] propose Mobile and Static sensor network Localization (MSL*) that extends MCL to the case where some or all nodes are static or mobile. On the other hand, Ref [3] provides Range-based Sequential Monte Carlo Localization method (Range-based SMCL), which combines range-free and range-based information to improve localization performance in mobile sensor networks. However, MCL methods are time-consuming as they need to keep sampling and filtering until enough samples are obtained to represent the posterior distribution of a mobile robot’s position[39, 16].

On the other hand, Extended Kalman Filter (EKF) approaches, [17, 18, 19, 20, 21, 22], provide suboptimal solutions by linearizing the measurements around the robot’s current position estimate, [21]. In particular, Ref. [17] provides a distributed EKF algorithm to estimate the location as well as the rotation of a robot by decomposing a single KF into a number of smaller communication filters. Improvements based on Covariance Intersection (CI) is provided in [18, 19, 20]. The problem of relative localization of two mobile robots in $\mathbb{R}^{2}$ is considered in [22]. By applying the observability rank condition for nonlinear systems, the authors in  [22] show that the underlying physical nonlinear system of Cooperative Localization (CL) in general has three unobservable degrees of freedom, corresponding to the global position and orientation, i.e., the system is observable only for more than three landmarks. When both robots move simultaneously, the error grows unbounded and the configuration is not observable.

Relevant to our work, are [34] and [35]; Ref. [34] addresses CL in 3D using combinations of range and bearing measurements in conjunction with ego-motion measurements. Ref. [35] on the other hand, proposes a decentralized method to perform mutual localization in multi-robot systems using relative measurements that do not include the identity of the measured robot.

In this paper, we provide a distributed algorithm to solve global localization in a network of mobile robots, i.e., we assume that there is no prior estimates of the initial locations. As a motivating scenario, consider a multi-robot network of ground/aerial vehicles with no central or local coordinator and with limited communications, whose task is to transport goods in an indoor facility, where GPS signals are not available. In order to perform a delivery task, each mobile robot has to know its own location first. In such settings, we are interested in developing distributed algorithms to track the robot locations such that the convergence is invariant to the initial position estimates. However, the implementation of distributed algorithms in mobile networks, where both robots and known beacon(s) are moving, is not straightforward due to the following challenges:

(i) no robot may be in proximity of any known beacon (hereinafter, referred to as a beacon for simplicity);

(ii) a robot may not be able to find nearby robots at all times to perform a distributed algorithm; and, (iii) the dynamic neighborhood at each robot results into a time-varying distributed algorithm, whose stability (convergence) analysis is non-trivial.

We address these challenges by providing an opportunistic update scenario, where a robot updates its location estimate in $\mathbb{R}^{m}$ only if it lies inside a convex hull of $m+1$ neighbors. Such neighbors are referred to as a triangulation set. Using this approach, we show that the robot location estimates are improved as the procedure continues and the algorithm is optimal, i.e., it tracks the true robot locations. In this context, the main contribution of this work is to develop a linear framework for localization that enables us to circumvent the challenges posed by the predominant nonlinear approaches to this problem. This linear framework is not to be interpreted as a linearization of an existing nonlinear algorithm. Instead, the nonlinearity from range to location is embedded in an alternate representation provided by the barycentric coordinates.

We abstract our localization algorithm as an LTV system, whose system matrices may be stochastic or sub-stochastic. These system matrices do not belong to a finite or a countable set rendering many of the existing results inapplicable. In addition, since they can be either stochastic or sub-stochastic, their product is no longer a group under multiplication. Thus, establishing the asymptotic behavior of the underlying LTV system is non-trivial. To address these challenges, we apply a novel method, [40, 41], to study LTV convergence by partitioning the entire set of system matrices into slices and relating the convergence rate to the slice lengths. In particular, we show that the algorithm converges if the slice lengths do not grow faster than a certain exponential rate. Next, since our localization scheme is based on the motion and the distance measurements, it is meaningful to evaluate the performance of the algorithm when these parameters are corrupted by noise. Hence, we study the impact of noise on the convergence of the algorithm and provide modifications to counter the undesirable effects of noise. Finally, we investigate what role does motion play in the behavior of the localization algorithm, and provide necessary conditions in terms of the number of nodes and the dimensions of motion required by our approach. As a consequence, the proposed method does not work, e.g., in the case of one beacon and less than three robots (in $\mathbb{R}^{2}$). Alternative approaches (such as SLAM and filtering-based localization methods) may be used in such settings, or in scenarios where a robot is isolated from the entire network.

The rest of this paper is organized as follows. In Section II, we formulate the problem. We propose our localization algorithm in Section III, followed by the convergence analysis in Section IV. We investigate the effects of noise in Section V, while in Section VI we discuss the relationship between the dimension of motion and the minimum number of beacons required. We discuss different aspects of the algorithm in Section VIII and present detailed simulation results in Section VII. Finally, Section IX concludes the paper.

II Preliminaries and problem formulation

Consider a network of $N$ robots with unknown locations in the index set, $\Omega$, and $M$ beacons with known locations in the index set, $\kappa$. Let $\Psi=\Omega\cup\kappa$ be the set of all nodes (robots and beacons), possibly mobile, located in $\mathbb{R}^{m},~{}m\geq 1$. We assume that the robots occupy a non-trivial configuration, i.e., all of them do not remain on a low-dimensional subspace in $\mathbb{R}^{m}$ effectively reducing the localization problem to $\mathbb{R}^{m-1}$. We denote the true location of the $i$-th node, $i\in\Psi$, at time $k$, with an $m$-dimensional row vector, $\mathbf{x}_{k}^{i\ast}\in\mathbb{R}^{m}$, where $k\geq 0$ is the discrete-time index. The problem is to find the locations of the mobile robots in the set $\Omega$, given any initialization of the underlying algorithm. In what follows, we explain our system model and describe a related convex hull inclusion test. We then provide a set of assumptions, which we use later in the design and analysis of our algorithm.

II-A System model

We express the motion as the deviation from the current to the next locations, i.e.,

[TABLE]

in which ${\widetilde{\mathbf{x}}_{k}^{i}}$ is the true motion vector at time $k$. We assume that robot $i$ measures a noisy version, ${\widehat{\mathbf{x}}_{k}^{i}}$, of this motion, e.g., by integrating the measurements gathered by an accelerometer:

[TABLE]

where $n^{i}_{k}$ represents the integration noise at time $k$. We define the true distance between any two nodes, $i$ and $j$, at the time of communication, $k$, as $\widetilde{d}^{ij}_{k}$. We assume that the motion is restricted in a bounded region in $\mathbb{R}^{m}$ and that the measured distance, $\widehat{d}^{ij}_{k}$, at robot $i$ includes noise, i.e.,

[TABLE]

where $r^{ij}_{k}$ is the noise in the distance measurement at time $k$.

II-B Inclusion test and Barycentric representation

As we will explain, to perform an update, each robot, $i$, must find $m+1$ neighbors, say $\{1,2,3\}$ in $\mathbb{R}^{2}$, such that it lies strictly in the interior of their convex hull, $\mathcal{C}(\cdot)$–a triangle formed by $\{1,2,3\}$. These neighbors, $\{1,2,3\}$, form what we refer to as a triangulation set, $\Theta_{i}(k)$, at time $k$. We note that the triangulation set, $\Theta_{i}(k)$, has a non-zero volume, i.e., $A_{\Theta_{i}(k)}>0$, where $A_{\Theta_{i}(k)}$ represents the $m$-dimensional volume, area in $\mathbb{R}^{2}$ or volume in $\mathbb{R}^{3}$, of $\mathcal{C}(\Theta_{i}(k))$, see Fig. 1, and can be computed by using the Cayley-Menger (CM) determinant, [42], see Appendix A. Clearly, $i\in\mathcal{C}(\Theta_{i}(k))$ and $|\Theta_{i}(k)|=m+1$. Thus to perform an update, a robot implements an inclusion test to verify if it lies inside an available convex hull. In any arbitrary dimension, $m$, a convex hull inclusion test can be described as follows:

[TABLE]

in which $``\setminus"$ denotes the set difference. Straightforward computation of Cayley-Menger determinants in $\mathbb{R}^{2}$ and $\mathbb{R}^{3}$, allows us to easily compute and take advantage of the barycentric coordinates in order to avoid nonlinearity in the localization process. CM determinant provides a simple equation that extends to any dimensions. It computes the area/volume of m-dimensional simplexes and reduces to, e.g., Heron’s formula for the area of a triangle in $\mathbb{R}^{2}$, and volume of a tetrahedron in $\mathbb{R}^{3}$. We note that $A_{\Theta_{i}(k)\cup\{i\}\setminus j}$ is the volume (area) of the set $\Theta_{i}(k)$ with node $i$ added and node $j$ removed. Furthermore, the test only requires pairwise distances among the nodes in $\{i\}\cup\mathcal{C}_{i}(\cdot)$, see Section VIII-B for the associated computation complexity.

With the help of triangulation sets, we represent the location of node $i$ at time $k$, ${\mathbf{x}_{k}^{i\ast}}$, as a linear-convex combination of the locations of the neighboring nodes in its triangulation set, see Fig. 1 (Left):

[TABLE]

where $a_{k}^{ij}$’s are the barycentric coordinates, defined as

[TABLE]

The above representation dates back to the early work by Lagrange and Möbius, [43]. Since we assume robot $i$ to lie strictly inside its triangulation set, all barycentric coordinates are strictly positive, i.e., $a_{k}^{ij}>0,\forall j\in\Theta_{i}(k)$.

II-C Assumptions

We now enlist our assumptions:

A0: Beacon locations $\mathbf{u}_{k}^{m\ast},m\in\kappa$, are known at all times.

A1: Each robot, $i\in\Omega$, knows a noisy version, $\widehat{\mathbf{x}}_{k}^{i}$, of its motion, $\widetilde{\mathbf{x}}_{k}^{i}$, at all times, see Eq. (2).

A2: Each robot, $i\in\Omega$, measures a noisy distance, $\widehat{d}^{ij}_{k}$, to every node, $j\in\Psi$, in its communication radius, $r$, at time $k$, see Eq. (3).

We assume that each robot is able to estimate its distances to the nearby nodes by using Received Signal Strength Indicator (RSSI), Time of Arrival (ToA), Time Distance of Arrival (TDoA), Direction of Arrival (DoA) [44], or camera-based methods, [45, 46]. Under the above assumptions, we are interested in finding the true locations of each robot without the presence of any central or local coordinator. In the following, Sections III-IV, we consider the ideal scenario when the motion and distance measurements are not effected by noise, i.e., $n_{k}^{i}=0$ and $r_{k}^{ij}=0$, in Eqs. (2) and (3). We then study the effects of noise in Section V, and provide modifications to the algorithm to address the noise on the motion as well as the distance measurements.

III Localization Algorithm

Consider a network of $|\Omega|=N$ mobile robots and $|\kappa|=M$ possibly mobile beacons in $m$-dimensional Euclidean space, ${\mathbb{R}}^{m}$. We define ${\mathcal{N}}_{i}(k)\subseteq\Psi$ as the set of neighbors of robot, $i\in\Omega$, at time $k$. We now describe our localization algorithm (main steps of the proposed localization algorithm are summarized in Algorithm 1): At the beginning, each robot starts with a random guess of its location estimate (Alg. 1 line 2). At each time $k>0$, we consider the following update scenarios for any arbitrary robot, $i\in\Omega$:

Case (i): If robot $i$ does not find at least $m+1$ neighbors, i.e., $0\leq|{\mathcal{N}}_{i}(k)|<m+1$, it does not update its current location estimate (Alg. 1 lines 7 and 8).

Case (ii): If robot $i$ finds at least $m+1$ neighbors, i.e., $|{\mathcal{N}}_{i}(k)|\geq m+1$, it performs the inclusion test, described in Section II-B, on all possible combinations of $m+1$ neighbors (Alg. 1 line 10). If robot $i$ then fails to find a triangulation set, it does not update (Alg. 1 line 12), otherwise it applies the following update111We consider two update scenarios for case (ii); A robot can update with respect to all possible triangulation sets at time $k$, or alternatively, it can only update once as soon as the inclusion test is passed for the first time at time $k$. (Alg. 1 line 14):

[TABLE]

where $\mathbf{x}^{i}_{k}$ is the location estimate of robot $i$ at time $k$, $\Theta_{i}(k)$ is the triangulation set at time $k$, {\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}a_{k}^{ij}>0} is the barycentric coordinates of node $i$ with respect to the node $j\in\Theta_{i}(k)$, and $\alpha_{k}$ is such that

[TABLE]

in which $\beta$ is a design parameter, see Section IV-C (v).

A network of $|\Omega|=5$ mobile robots and $|\kappa|=2$ beacons with fixed locations in $\mathbb{R}^{2}$ is illustrated in Fig. 2. To clarify the updating scenarios let us just consider robot $1$, shown as a blue filled circle. At time $k_{j}$, Fig. 2 (a), robot $1$ has only one neighbor and does not update, as per case (i). Although at time $k_{l}$, Fig. 2 (b), robot $1$ finds $m+1=3$ nodes (two robots and one beacon) within its communication radius, it fails to update because the neighbors do not satisfy the inclusion test. In this case, $|\mathcal{N}_{1}(k_{l})|=3$ and $\Theta_{1}(k_{l})=\emptyset$. However, robot $1$ updates its location estimates at time $k_{m}$ and $k_{n}$, where in the former it finds two different triangulation sets: $\{2,3,4\}$ and $\{2,4,5\}$, and disregards $\{3,4,5\}$ and $\{2,3,5\}$.

At any given time, $k$, note the difference between $\mathcal{N}_{i}(k)=\emptyset$ and $\Theta_{i}(k)=\emptyset$; while the former implies that node $i$ has no neighbors, as in case (i), the latter implies that robot $i$ cannot find any triangulation set among the neighbors.222Note that our algorithm can be generalized to higher dimensions, e.g., $3$-D, where triangles and areas become tetrahedrons and volumes, respectively. For example at time $k_{j}$, Fig. 2 (a), we have $\mathcal{N}_{i}(k)=\{2\}$ and $\Theta_{i}(k)=\emptyset$, while at time $k_{n}$, Fig. 2 (d), we have $\mathcal{N}_{i}(k)=\Theta_{i}(k)=\{3,4,5\}$.

By separating the barycentric coordinates over the robots in $\Omega$ and the beacons in $\kappa$, we can split Eq. (8) as

[TABLE]

where $\mathbf{u}_{k}^{m}$ is the true location of $m$-th beacon at time $k$ and

[TABLE]

It can be inferred from Eqs. (11) and (12) that at time $k$ the self-weight assigned to the $i$-th robot’s estimate, $p_{k}^{ii}$, is always lower bounded, i.e.,

[TABLE]

Although it is possible that more than one robot finds a triangulation set and update at time $k$, in the remaining of the paper we assume, without loss of generality, that at most one robot updates at each iteration. We can now write the above algorithm in matrix form as follows:

[TABLE]

in which ${\bf{x}}_{{k}}\in\mathbb{R}^{N\times m}$ is the vector of robot coordinates at time $k$, ${\bf{u}}_{{k}}\in\mathbb{R}^{M\times m}$ is the vector of beacon coordinates at time $k$, and $\widetilde{{\bf{x}}}_{{k+1}}\in\mathbb{R}^{N\times m}$ is the change in the location of robots at the beginning of the $k$-th iteration according to Eq. (1). Note that ${\bf{P}}_{{k}}$, and ${\bf{B}}_{{k}}$, the system matrix and the input matrix of the above LTV system, contain (weighted) barycentric coordinates at time $k$, with respect to the robots with unknown locations, and beacons, respectively.

Since true information is only injected into the network by the beacons, we must set a lower bound on the weights assigned to the beacon states. Otherwise, the beacons may be assigned a weight that goes to zero over time, i.e., the beacons eventually are excluded from the network. To make sure that beacons remain in the network, we naturally make the following assumption.

A3: Guaranteed beacon contribution. When a beacon is involved in an update, i.e., for any $b_{k}^{im}\neq 0$, we impose that

[TABLE]

where $\alpha$ is the minimum beacon contribution, see Section IV-C (vi).

The above assumption implies that if there is a beacon in the triangulation set, a certain amount of information is always contributed by the beacon. In other words, the update occurs only if the robot lies in an appropriate location inside the convex hull. By Assumption A3, the weight assigned to the beacon should be at least $\alpha$. This weight comes from the barycentric coefficient corresponding to the beacon, i.e.,

[TABLE]

where $i$ and $m$ indices represent the updating robot and a beacon, respectively. To clarify Assumption A3, we consider Fig. 3 where the updating robot lies inside the convex hull of three nodes, including two other robots and one beacon; let $\alpha=0.25$. Fig. 3 (Left) shows a situation, where the robot updates, and the beacon provides the exact minimum contribution, i.e., the area of the shaded triangle is one fourth of the area of the convex hull triangle.

On the other hand, if the robot lies within the shaded trapezoid illustrated in Fig. 3 (Middle), Eq. (17) becomes less than $\alpha$, and Eq. (16) does not hold. No update occurs in this case since the beacon does not provide enough valuable information for the robot. The upper side of the shaded trapezoid is the threshold boundary, such that if the robot lies on this line, an update occurs and the beacon provides the exact minimum contribution. Thus, our measure to evaluate the contribution of a beacon in an update, is whether or not the updating robot lies inside a specific area, and not the distance between the robot and the beacon. If the updating robot stays inside the same convex hull, and moves closer to the beacon, as shown in Fig. 3 (Right), the robot performs an update and the weight assigned to the beacon exceeds $\alpha$. Note that in this example we set $\alpha=0.25$ only for the ease of demonstration. The localization algorithm works for any arbitrary small (non-zero) value of $\alpha$, which corresponds to higher possibilities of a successful update when a beacon is involved.

With the lower bounds on the self-weights according to Eq. (14), and the lower bounds on the weights assigned to the beacons according to Eq. (16), the matrix of barycentric coordinates with respect to robots with unknown locations, i.e., the system matrix at time $k$, ${\bf{P}}_{k}$, is either

1. (i)

identity, when no robot updates at time333Note that this case also characterizes the case, where a robot fails to update due to a (temporary) communication loss/drop. $k$; or, 2. (ii)

identity except a stochastic $i$-th row, when a robot, say $i$, updates but finds no beacon in its triangulation set at time $k$, i.e.,$\Theta_{i}(k)\neq\emptyset$ and $\Theta_{i}(k)\cap\kappa=\emptyset$; or, 3. (iii)

identity except a sub-stochastic $i$-th row, when the updating robot, say $i$, has at least one beacon in its triangulation set at time $k$, i.e., $\Theta_{i}(k)\cap\kappa\neq\emptyset$.

Note that when a robot, say $i$ at time $k$, performs an update only with robots in the set $\Omega$ (case (ii) above), the $i$-th row of the system matrix, $\mathbf{P}_{k}$, contains all of the non-zero weights (corresponding to the barycentric coordinates) and $\mathbf{B}_{k}$ is a zero matrix. Each row of $\mathbf{P}_{k}$ is identity except the $i$th row whose sums is one due to convexity. However, if a beacon in the set $\kappa$ is involved in an update (case (iii) above), the $i$-th row of the input matrix, $\mathbf{B}_{k}$, contains a non-zero weight assigned to that beacon. Therefore the weights assigned to the other neighboring robots sum to a value strictly less than one; in other words, $\mathbf{P}_{k}$ is identity except the $i$-th row whose sum is strictly less than one making $\mathbf{P}_{k}$ sub-stochastic.

In the next section, we investigate the asymptotic stability of LTV systems where the system matrices follow the cases described above, and provide sufficient conditions for the iterative localization algorithm, Eq. (15), to converge to the true robot locations. Before we proceed, let us make the following definitions to clarify what we mean by stochasticity and sub-stochasticity throughout this paper:

Definition 1**.**

A non-negative, stochastic matrix is such that all of its row sums are one. A non-negative, sub-stochastic matrix is such that it has at least one row that sums to strictly less than one and every other row sums to at most one.

IV Convergence analysis

We start this section by providing a related result on the convergence of an infinite product of stochastic/sub-stochastic matrices, [40]. We will then adapt these results to investigate the convergence of our localization algorithm represented by Eq. (15).

IV-A Asymptotic behavior

Consider the following LTV system:

[TABLE]

in which the time-varying system matrix, $\mathbf{P}_{k}\triangleq\{p_{k}^{ij}\}$, is non-negative, random, and represents at most one state update, say at the $i$-th row, for any $k$. In other words, $\mathbf{P}_{k}$’s can randomly switch between stochastic, and sub-stochastic matrices. In addition, assume the following on the update at time $k$:

B0: If $\sum_{j}p_{k}^{ij}=1$, i.e., the updating row, $i$, in $\mathbf{P}_{k}$ is stochastic, then the $i$-th self-weight is lower-bounded as

[TABLE]

B1: If $\sum_{j}p_{k}^{ij}<1$, i.e., the updating row, $i$, in $\mathbf{P}_{k}$ is sub-stochastic, then it is also upper-bounded as

[TABLE]

To study the asymptotic behavior of an LTV system with such system matrices, in [40], we introduce the notion of a slice, $M_{j}$, as the smallest product of consecutive system matrices, such that: (i) the infinity norm of each slice is less than one, i.e., $\|M_{t}\|_{\infty}<1,\forall t$; and, (ii) the entire sequence of system matrices is covered by non-overlapping slices, i.e.,

[TABLE]

Each slice initiates with a sub-stochastic system matrix, and terminates after all row sums become less than one, i.e., each row in the slice becomes sub-stochastic. The length of the $j$-th slice, ${|{M_{j}}|}$, is defined as the number of matrices forming the slice, and the upper bound on the infinity norm of a slice is further related to its length, [41], as

[TABLE]

The following theorem characterizes the asymptotic behavior of the LTV system represented by Eq. (18):

Theorem 1**.**

*With Assumptions B0-B1, the LTV system, $\mathbf{x}_{k+1}=\mathbf{P}_{k}\mathbf{x}_{k}$, converges to zero , i.e., $\lim_{k\rightarrow\infty}\mathbf{x}_{k}=\mathbf{0}_{N}$, if any one of the following is true: *

1. (i)

Each slice has a bounded length, i.e.,

[TABLE] 2. (ii)

There exist an infinite subset, $J_{1}$, of slices such that**

[TABLE] 3. (iii)

For every $i\in\mathbb{N}$, there exists an infinite subset, $J_{2}$, of slices such that**

[TABLE]

for some $\gamma_{1}\in[0,1]$, $\gamma_{2}>0$. For any other slice, $M_{j},~{}j\notin J_{2}$ we have $|M_{j}|<\infty$.

We now briefly discuss the intuition behind the above theorem. According to Eq. (21), we can use the product of slices instead of the product of system matrices, and study the following dynamics rather than Eq. (18): ${\bf{y}}({t+1})={\bf{M}}_{t}{\bf{y}}({t})$. By taking the infinity norm of the both sides of this equation and using the sub-multiplicative norm property we get

[TABLE]

The first two cases are trivial, because if all or an infinite subset of slices have a bounded length, then according to Eq. (22) each slice will have a bounded (and subunit) infinity norm, whose infinite product is zero. Thus, from Eq. (27) we can infer that ${\|{\bf{y}}(t+1)\|}_{\infty}=0$, which in turn leads to $\lim_{k\rightarrow\infty}\mathbf{x}_{k}=\mathbf{0}_{N}$, and completes the proof. However, we can show that a strict upper bound on the lengths of all or an infinite subset of slices is not necessary. In fact, all we require is an infinite subset of slices whose (unbounded) lengths grow slower than the exponential rate provided in Eq. (26). Later in this section, we explain how the slice lengths can be interpreted in the context of inter-network communications. The reader is referred to our prior works, [40, 41], for a detailed proof of Theorem 1. In what follows, we use the results of this theorem to study the convergence of Eq. (15).

IV-B Convergence of the localization algorithm

We start this section with the following lemma:

Lemma 1**.**

Under Assumptions A0-A3 and no noise, the product of system matrices, $\mathbf{P}_{k}$’s, in the LTV system represented by Eq. (15), converges to zero if any one of the three conditions in Theorem 1 holds.

Proof.

We need to show that B0 and B1 can be inferred from Assumptions A0-A3. First note that Eq. (14) results in B0 if we set $\beta_{1}=\beta$. On the other hand, if Assumption A3 holds, i.e., if there exist at least one beacon in the triangulation set, we can write

[TABLE]

in which we used the fact that barycentric coordinates sum to one. Note that the second term on the right hand side (RHS) of Eq. (28) gives the sum of all weights assigned to the beacon(s) in the triangulation set. This term is minimized, and hence the RHS of Eq. (28) is maximized, when there is only one beacon among the neighbors of the updating robot, and the minimum weight, $\alpha$, is assigned to this beacon according to Eq. (16). In this case we can write Eq. (28) as

[TABLE]

which provides an upper bound on the $i$-th row sum of ${\bf{P}}_{{k}}$, and results in Assumption B1 if we choose $\beta_{2}=1-\alpha$. With Assumptions B0-B1 satisfied, the rest of the proof follows from that of Theorem 1. ∎

We now provide our main result in the following theorem.

Theorem 2**.**

Consider a network of $M$ (possibly mobile) beacons and $N$ mobile robots moving in a finite and bounded region. Then, under Assumptions A0-A3 and no noise, for any (random or deterministic) motion that satisfies one of the conditions in Theorem 1, the solution of Eq. (15) asymptotically converges to the true robot locations.

Proof.

The motion of the robots in a finite, bounded region results in the following LTV system:

[TABLE]

in which the system matrices, ${\bf{P}}_{{k}}$’s switch between stochastic and sub-stochastic matrices. Therefore, under Assumptions A0-A3 and according to Lemma 1, if any of the conditions in Theorem 1 is satisfied, we have

[TABLE]

On the other hand, the true robot locations are given by:

[TABLE]

To justify the above equation, note that if robot $i$ lies inside the convex hull of $m+1$ neighbors at time $k$, its true location can be expressed as a convex combination of the true locations of its neighbors. Note that when the inclusion test is not passed at any robot at time $k$, the system matrix, ${\bf{P}}_{{k}}$, and the input matrix, ${\bf{B}}_{{k}}$, become identity and zero matrices, respectively, and the location estimate, ${\bf{x}}_{{k+1}}$, and the true locations, ${\bf{x}}^{*}_{{k+1}}$ in Eqs. (30) and (32) update according to the motion vector. By subtracting Eq. (30) from Eq. (32) the error dynamics can be obtained as follows

[TABLE]

which converges to zero from Lemma 1 and Eq. (31), and the proof is complete. ∎

IV-C Discussion

In what follows, we shed some light on Theorem 1, and elaborate on the choice of design parameters, $\alpha$ and $\alpha_{k}$.

(i) As mentioned earlier in Section IV-A, each slice is initiated with a sub-stochastic update, which only occurs when a beacon is among the neighbors of the updating robot. In other words, to initiate a slice a robot with unknown location must have a beacon in its triangulation set. Thus, the first system matrix of each slice has one sub-stochastic, and $N-1$ stochastic rows. We call a robot informed after it updates with a beacon. Therefore, at the beginning of a slice there are one informed and $N-1$ uninformed robots in the network.

(ii) A slice is terminated after all rows of a slice are sub-stochastic, i.e., all robots are informed. For this to happen, each robot has to either update with a beacon directly, or indirectly, i.e., update with a robot who updated with a beacon in the same slice. By indirectly, we mean that a robot receives information not from a beacon but by any other informed robot in the network. Completion of each slice corresponds to the propagation of location information from beacon(s) to all robots in the network, and the location estimates are refined each and every time a slice is completed. Hence, slice representation is closely related to the information flow from the beacons to every other robot in the network.

(iii) Theorem 1 provides the conditions on the rate, at which such information should propagate for convergence to the true robot locations. The first two cases in Theorem 1 require all or an infinite subset of slices to have bounded lengths. In other words, the information from the beacons has to reach all robots within $L$ or $L_{1}$ iterations, infinitely often. These conditions are relaxed in the third case; Eq. (26) implies that all we need is an infinite subset of slices whose lengths grow slower than a certain exponential rate. The condition on motion described in Theorem 2 can be translated to a condition on the information dissemination from Theorem 1. Thus, any motion that guarantees this information dissemination suffices.

(iv) The non-zero self-weights, $\alpha_{k}$’s, assigned to the previous state of the updating robot guarantees that an informed robot does not become uninformed again within the same slice, e.g., by performing an update with a set of uninformed robots before the slice is complete.

(v) Although the proposed algorithm converges for any value of $0<\alpha<1$, the convergence rate of the algorithm is affected by the choice of $\alpha$: By choosing $\alpha$ arbitrarily close to $1$, a robot has to get arbitrarily close to a beacon in order to perform an update with respect to that beacon (see Fig. 3), which in turn corresponds to arbitrary large number of iterations for the termination of each slice. On the other hand, setting $\alpha$ arbitrarily close to zero makes slice norms arbitrarily close to $1$. A proper choice can be made by considering the motion model, the communication protocol, and the number of available beacons in the network.

V Localization under imperfect measurements

The noise on the motion and distance measurements degrades the performance of the localization algorithm, as expected, and in certain cases the location error is larger than the region of motion; this is shown experimentally in Section VII. In what follows, we provide the modifications, M1-M3, to the proposed algorithm to counter the undesirable effects of noise in case of motion, ${\widehat{\mathbf{x}}_{k}^{i}}$, and on the distance measurements, $\widehat{d}^{ij}_{k}$.

V-A Discard unreliable Cayley-Menger determinants

If a robot is located close to the boundaries of a convex hull, the noise on distance measurements may affect the inclusion test results. To get meaningful values for the areas and volumes in $\mathbb{R}^{2}$ and $\mathbb{R}^{3}$ according to Eq. (40), the corresponding Cayley-Menger determinants computed with perfect distance measurements must obey a certain sign in order for the square root to convey a meaningful result (negative in $\mathbb{R}^{2}$ and positive in $\mathbb{R}^{3}$). Therefore, we suggest the first modification to the algorithm as follows:

M1: A robot does not perform an inclusion test if the corresponding Cayley-Menger determinant is positive in $\mathbb{R}^{2}$, or negative in $\mathbb{R}^{3}$.

V-B Inclusion test error

Even in the case of perfect distance measurements, the inclusion test results may not be accurate due to the noise on the motion, which in turn corresponds to imperfect location updates at each and every iteration. To address this issue, we propose the following modification to the algorithm:

M2: Suppose the inclusion test is passed at time $k$ by a triangulation set, $\Theta_{i}(k)$. Robot $i$ performs an update only if

[TABLE]

where $\epsilon_{k}^{i}$ is the relative inclusion test error at time $k$ for robot $i$, and $\epsilon$ is a design parameter. The optimal choice of $\epsilon$ depends on the number of robots and beacons in the network and the statistics of the noise.

V-C Convexity

Finally, in order to guarantee the convexity in the updates in the presence of noise, we consider the following modification:

M3: Suppose the inclusion test is passed by a triangulation set, $\Theta_{i}(k)=\{j,~{}l,~{}m\}$, and Eq. (34) holds at time $k$. Robot $i$ first computes, $a_{k}^{{i}j}$, $a_{k}^{{i}l}$, and $a_{k}^{{i}m}$, according to Eq. (7) and using the noisy distance measurement. It then normalizes the weights assigned to each neighbor in order to preserve the convexity of the update. For example, the (normalized) weight assigned to robot $j$ is computed as $a_{k}^{{i}j}/\sum_{n\in{\Theta_{i}(k)}}a_{k}^{{i}n}$.

In Section VII, we show that the above modifications to the algorithm improve the results significantly in presence of noise, and result in a bounded error in location estimates.

VI How many beacons are necessary?

In this section, we investigate the minimal number of beacons required for localizing an arbitrary number of robots using the algorithm described in Section III. We denote the subspace of motion at robot, $i\in\Omega$, and beacon, $j\in\kappa$, by $\mathcal{M}_{i}$ and $\mathcal{U}_{j}$, respectively. Let us clarify this notation with a simple example: Suppose robot $1$ is moving along a vertical line (regardless of the direction); this line forms $\mathcal{M}_{1}$, and $\dim\mathcal{M}_{1}=1$. Note that $\mathcal{M}_{i}$ or $\mathcal{U}_{j}$ includes all possible locations that the $i$-th robot or the $j$-th beacon occupies throughout the localization process, i.e., discrete times $k=1,~{}2,\ldots$. Now consider another robot, $2$, which is moving along a vertical line parallel to $\mathcal{M}_{1}$ (regardless of the direction); in this case we will have $\dim\cup_{i=1,2}\mathcal{M}_{i}=\dim\mathcal{M}_{2}=1$. However, if the two lines are linearly independent, they span $\mathbb{R}^{2}$, and we will have $\dim\cup_{i=1,2}\mathcal{M}_{i}=2$.

Assuming $\mathbb{R}^{m}$, we show that the motion of the robots and beacons in $l\leq m$ dimensions allows us to reduce the number of beacons from $m+1$ by $l$. Note that the traditional trilateration scheme requires at least $3$ nodes with known locations in $\mathbb{R}^{2}$. Therefore, assuming $m+1$ beacons in $\mathbb{R}^{m}$ has been standard in all multilateration-based localization algorithms in the literature, see e.g., [27, 47]. In the following theorem, we develop necessary conditions for the proposed algorithm to track the true location of robots. Since we can provide robots with up to $m$ degrees of freedom in their motion in $\mathbb{R}^{m}$, we then show that our localization algorithm works in the presence of only one ($m+1-m$) beacon.

Theorem 3**.**

For the LTV dynamics in Eq. (15) to track the true robot locations, following a non-trivial configuration in $\mathbb{R}^{m},m=2$, the following conditions must be satisfied

[TABLE]

Proof.

Eq. (35) is trivial, because when there is no beacon in the network, all system matrices become stochastic, and the error dynamics in Eq. (33) do not converge to zero. This is equivalent to Case (ii) in Section III, where the update is always in terms of robots, i.e., $\mathbf{B}_{k}$ is always zero in Eq. (15). Eq. (36) stems from the fact that each robot requires at least $m+1$ neighbors to perform an update in $\mathbb{R}^{m}$, assuming non-trivial configurations. In order to prove the necessary condition in Eq. (37), let us first consider the $2$-dimensional Euclidean space. Suppose on the contrary that Eq. (37) does not hold. In this case, there exist three possible scenarios as follows:

1. (i)

$|\kappa|=0,\mbox{and}\dim\underset{i\in\Omega}{\cup}\mathcal{M}_{i}=2$; or, 2. (ii)

$|\kappa|=1,\mbox{and}~{}\dim\underset{i\in\Omega}{\cup}\mathcal{M}_{i}+\dim\underset{j\in\kappa}{\cup}\mathcal{U}_{j}=0$; or, 3. (iii)

$|\kappa|=1,\mbox{and}~{}\dim\underset{i\in\Omega}{\cup}\mathcal{M}_{i}+\dim\underset{j\in\kappa}{\cup}\mathcal{U}_{j}=1$; or, 4. (iv)

$|\kappa|=2,\mbox{and}~{}\dim\underset{i\in\Omega}{\cup}\mathcal{M}_{i}+\dim\underset{j\in\kappa}{\cup}\mathcal{U}_{j}=0$.

Localization is not possible in case (i), as it violates Eq. (35).

Case (ii)–with $1$ beacon, $3$ robots and no motion, is where at most one robot may be able to lie inside the convex hull of the remaining three; making the triangulation of the second robot impossible. Same argument can be applied for Case (iv)–with $2$ beacons, $2$ robots and no motion. Mathematically, these two cases mean that a slice is never completed because the rows of $\mathbf{P}_{k}$ in Eq. (15) corresponding to the non-updating robots will always be the corresponding rows of identity. Thus, all we need to show is that localization is not possible in case (iii), where there is one beacon in the network, and the dimension of the motion in all nodes is one, i.e., all nodes can only move along parallel lines in $\mathbb{R}^{2}$.

We start with the best possible scenario, which initially allows one robot to triangulate.

We then choose an arbitrary direction in which all nodes are allowed to move. As shown in Fig. 4 (Left), the beacon and the three robots can move along $\ell$, $\ell_{1}$, $\ell_{2}$, and $\ell_{3}$, respectively. If this direction is chosen such that the beacon’s line of motion lies within the strip between $\ell_{1}$ and $\ell_{2}$ (the lines of motion for the two robots that are not initially able to triangulate), no other robot will be able to triangulate. That is because robots $1$ and $2$ have to move inside the shaded regions illustrated in Fig. 4 (Right) in order to (possibly) triangulate. Clearly, any vector that provides such motion is linearly independent of the robots’ direction of motion.

We now consider the case where the direction of motion is such that the beacon lies outside the aforementioned strip. This scenario is illustrated in Fig. 5.

At time $k_{2}$, robot $3$ lies inside the convex hull of the beacon and the other two robots, and is able to triangulate, see Fig. 5 (Left). As shown in Fig. 5 (Right), the robots can move along the specified parallel lines such that at time $k_{3}$ robot $1$ lies inside the triangle of the beacon and robots $2$ and $3$. However, the farthest robot from the beacon (robot $2$ in this case) will never be able to triangulate, unless it moves inside the shaded region shown Fig. 5 (Right), i.e., the strip between the lines of motion corresponding to robot $3$ and the beacon. This in turn requires a motion vector that is linearly independent of the robots’ specified direction of motion, i.e., $\dim\underset{i\in\Omega}{\cup}\mathcal{M}_{i}$ greater than $1$. At time $k_{2}$, Fig. 5 (Left), robot $3$ lies inside the convex hull of the beacon and the other two robots, and is able to triangulate. Due to the motion, Fig. 5 (Right), the robots move along the specified parallel lines such that at time $k_{3}$, robot $1$ lies inside the triangle of the beacon and robots $2$ and $3$. However, the farthest robot from the beacon (robot $2$ in this case) will never be able to triangulate, unless it moves inside the shaded region shown in Fig. 5 (Right), i.e., the strip between the lines of motion corresponding to robot $3$ and the beacon. This in turn requires a motion vector that is linearly independent of the robots’ specified direction of motion, i.e., $\dim\underset{i\in\Omega}{\cup}\mathcal{M}_{i}$ greater than $1$. Mathematically, row $2$ of $\mathbf{P}_{k}$ in Eq. (15) never changes and a slice never completes. ∎

Although Theorem 3 is applicable to $m=3$, a formal proof is beyond the scope of this paper and will be provided elsewhere. In the sequel, we only give a sketch of the proof to show that localization is not possible with one beacon and two dimensional motion in $\mathbb{R}^{3}$. This scenario is shown in Fig. 6, where robot $4$ initially lies inside the convex hull (a tetrahedron in $\mathbb{R}^{3}$) of the other nodes.

The beacon and robots $1-4$ can move on the planes $\mathcal{P}$, and $\mathcal{P}_{1}-\mathcal{P}_{4}$ as shown in Fig. 6. Clearly, the farthest robot from the beacon, robot $2$ in this case, is not able to find a triangulation set, unless it moves into the space between $\mathcal{P}_{3}$ and $\mathcal{P}$, which in turn requires a motion vector which is not in the span of any pair of linearly independent vectors in the motion planes, and thus makes localization of robot $2$ impossible.

Theorem 3 provides necessary conditions for our localization algorithm in terms of the minimal number of beacons and the dimension of motion in the network.

In what follows, we show how adding a new beacon or a new dimension in the motion make the localization algorithm work in $\mathbb{R}^{2}$. We start with case (ii) in Theorem 3 and add one more beacon, i.e., $|\kappa|=2,\mbox{and}~{}\dim\underset{i\in\Omega}{\cup}\mathcal{M}_{i}+\dim\underset{j\in\kappa}{\cup}\mathcal{U}_{j}=1$. This scenario is illustrated in Fig. 7. Without loss of generality, we assume that only one robot, $1$, is moving from time $k_{1}$ to $k_{2}$. To triangulate, this robot has to lie inside the convex hull of the two beacons and robot $2$, the shaded region in Fig. 7 (a). To this aim, one possible motion vector is shown in Fig. 7 (b) with a solid vector. Note that in order to let robot $2$ triangulate, robot $1$ can move via the same motion vector (and in the opposite direction) to lie inside the shaded region shown in Fig. 7.

We now revisit case (ii) in Theorem 3 and add one more motion dimension, i.e., $|\kappa|=1,\mbox{and}~{}\dim\underset{i\in\Omega}{\cup}\mathcal{M}_{i}+\dim\underset{j\in\kappa}{\cup}\mathcal{U}_{j}=2$.

A possible motion trajectory of one robot, $1$, which leads to the triangulation of all robots in a network of size $4$ is shown in Fig. 8, assuming that the beacon and all the robots, except, robot $3$, are static. We denote the motion vectors with ${\bf{v}_{1}}$ and ${\bf{v}_{2}}$, which are orthogonal in this case. The shaded regions in Fig. 8 (a), (b), and (c) represent the locations, where robots $4$, $1$, and $3$ can triangulate. For example, Robot $2$ lies inside a triangle of Robots $3$, Robot $1$, and the beacon, only if Robot $3$ moves to the shaded region in Fig. 8 (a) between time $k_{1}$ and $k_{2}$. Note that Fig. 8 illustrates a scenario, in which one slice is completed, i.e., each robot receives information from the beacon at least once. Clearly, to achieve the convergence to the true locations, such motion has to be repeated infinitely often such that one of the conditions in Theorem 1 is satisfied.

VII Simulations

In this section, we provide the simulation results to illustrate the proposed localization algorithm in $\mathbb{R}^{2}$. Let the true location of the $i$-th robot, $i\in\Omega$ in $\mathbb{R}^{2}$, be decomposed as ${\mathbf{x}_{k}^{i\ast}}=[{x_{k}^{i\ast}}~{}{y_{k}^{i\ast}}]$. We consider the following random motion model for robot $i\in\Omega$:

[TABLE]

in which $d^{i}_{k+1}$ and $\theta^{i}_{k+1}$ denote the distance and angle traveled by robot $i$, between time $k$ and $k+1$, and are random. We choose, $d^{i}_{k}$ and $\theta^{i}_{k}$, to have uniform distribution over the intervals of $[0~{}d_{\max}]$ and $[0~{}2\pi]$, respectively, such that each robot, $i\in\Omega$, does not leave the bounded region of interest. Also note that the random motion is assumed to be statistically independent among the robots. In what follows, we first provide the simulation results in noiseless scenarios, and then study the effects of noise and the proposed modifications in the convergence of the algorithm.

VII-A Localization in noiseless scenarios

We first consider a network of $5$ mobile robots with unknown locations, and only one beacon, whose location is fixed and perfectly known at all time.

In the beginning, all robots are randomly deployed inside a $20\mbox{m}\times 20\mbox{m}$ square.

We fix the location of the beacon in the center of the region. For all simulations in this section, we set the communication radius to $r=2$m and $d_{\max}=5\mbox{m}$. All robots are initially assigned with a random estimate of their initial locations. If a robot finds at least $m+1=3$ neighbors, it performs the inclusion test, as described in Section II-B. If the opportunity occurs, i.e., if a robot finds a triangulation set among the neighboring nodes, it updates its location estimate according to Eq. (12). The robot does not perform any update otherwise. Throughout this section, we consider the scenario, where a robot performs multiple updates at each iteration with respect to all possible triangulation sets found in that iteration. To ensure that the updating robot retains the valuable information it may have received from the beacon, and to guarantee a minimum contribution by the beacon when it is involved in an update, we set $\alpha_{k}=\beta=0.01$ and $\alpha=0.01$, respectively.

Fig. 9 (Left) shows the motion model that we choose as random according to Eq. (38), i.e., the trajectories of $N=5$ mobile robots for the first $20$ iterations. We choose the second norm of the error vector, $\mathbf{e}_{k}^{i}$, to characterize the convergence, i.e.,

[TABLE]

The algorithm converges to the true robot locations as ${\|\mathbf{e}_{k}^{i}\|_{2}}\rightarrow 0,\forall~{}i$. The convergence of the algorithm in this case is illustrated in Fig. 9 (Right) for one simulation. In Fig. 10 (Left), we provide the convergence results for the networks with one beacon and $5$, $10$, $20$ and $100$ robots with unknown locations. Each curve indicates the average over $n=20$ Monte Carlo simulations.

Due to the opportunistic nature of the algorithm, as the number of robots increases each robot is more likely to find triangulation sets and perform successful updates. Therefore, the algorithm converges faster as the number of robots increases. In Fig. 10 (Right), we study the effect of the self-weights on the convergence rate of the algorithm in networks with one beacon and $5$, $10$ and $20$ robots as we change $\alpha_{k}=\beta$ from $0.01$ to $0.5$. It can be seen that the convergence rate decreases as $\alpha_{k}$ increases.

In Fig. 11 (Left), we show the percentage of the iterations a robot finds different number of neighbors in networks with one beacon and $5$, $10$ and $20$ robots. For example, a robot finds no neighbors during $40\%$ of the iterations in a network of $5$ robots. Fig. 11 (Right) shows the ratio of the total number of updates to the number of iterations. We fix the number of iterations to $3000$, and take the average over $n=20$ Monte Carlo simulations. On average, a robot in networks with one beacon and $5$, $10$, and $20$ robots, performs $31$, $171$, and $422$ updates, respectively.

We examine the effect of highly adverse initial conditions in Fig. 12 (Left) for a network with one beacon and $100$ robots. Despite very large initial errors, which is $30$ times larger than the dimension of the region, the algorithm converges in less than $500$ iterations in this simulation. This convergence rate is slower than the average convergence rate illustrated in Fig. 10 (Left) for networks with $100$ robots primarily because of the large initial error.

VII-B Localization in the presence of noise

We use two different models to examine the effects of noise on the proposed localization algorithm;

First, we assume that the noise on odometry measurements, i.e., the distance and angle that robot $i$ travels at time $k$, are Gaussian with zero mean and the following variances

[TABLE]

where $D^{i}_{k}$ represents the total distance that robot $i$ has traveled up to time $k$. We also assume that the noise on the distance measurement (to a neighboring robot) at time $k$ is normal with zero mean and the variance of ${\sigma_{r}^{i}}^{2}={K_{r}}^{2}k$. Therefore the variances of the odometry measurements are proportional to the total distance a robot has traveled, and the variance on the distance measurements (to the neighboring robots) increases with time. Such assumptions are common in the relevant literature, e.g., [22, 48, 49, 50].

As shown in Fig. 12 (Left) for a network with one beacon and $100$ robots, setting ${K_{d}}=K_{\theta}=K_{r}=5*10^{-3}$ leads to an unbounded error, which is due to incorrect inclusion test results and the continuous location drifts because of the noise on the distance measurements and the noise on motion, respectively. However, by modifying the algorithm according to Section V, it can be seen in Fig. 13 that the localization error is bounded by the communication radius. In the simulations with noise we choose $\epsilon=20\%$, i.e., a robot performs an update only if the relative inclusion test error, corresponding to the candidate triangulation set is less than $20\%$.

We evaluate the performance of the algorithm on a different noise model, where at each and every iteration the amount of noise on odometry and distance measurements are proportional to the measurements. In Fig. 14, we show the simulation results when the amount of noise at each iteration is up to $\pm 5\%$ of the measurements. All the simulations in the presence of noise are averaged over $n=20$ Monte Carlo simulations.

VII-C Performance evaluation

We now evaluate the performance of our algorithm in contrast with some well-known localization methods; MCL [5], MSL* [6], and Range-based SMCL [3].

Please refer back to Section I for a brief description of these methods. In Fig. 15, we compare the localization error in the Convex Hull (CH) algorithm with MCL, MSL*, and Range-based SMC. As shown in Table I, in these algorithms node density, $n_{d}$, and beacon (seed) density, $n_{s}$, denote as the average number of nodes and beacons in the neighborhood of an agent, respectively. Total number of agents and beacons can therefore be determined by knowing these densities as well as the area of the region. We consider $N=100$ robots and $M=10$ beacons to remain consistent with the setup in [5, 6, 3], and use the same metric, i.e., the location error as a percentage of the communication range. Each data point in Fig. 15 is computed by averaging the results of $20$ simulation experiments. We keep the other parameters the same as described earlier in this section. With high measurement noise, i.e., in the presence of $10\%$ noise on the range measurements and $1\%$ noise on the motion, our algorithm outperforms MCL, MSL* and Range-based SMCL after $10$ iterations. Clearly, the localization error in our algorithm decreases as the amount of noise decreases, and our algorithm converges to the exact robot locations in the absence of noise. Table I summarizes the performance of the proposed CH algorithm in comparison to the above methods.

VIII Remarks

VIII-A Convergence rate

In some particular applications, e.g., search and rescue in hazardous environments, the robots need not only to find their locations, but also to finish the localization process successfully in a finite time. The convergence rate of the algorithm becomes crucial in such applications. Since the algorithm is asymptotic, one can design appropriate termination criteria that are application-dependent. For example, one such criterion can be designed according to the number of iterations typically needed given the size, mobility, models, and noise parameters, as evident from the simulation figures in Section VII. We do note that the convergence improves dramatically as the size of the network increases.

VIII-B Computational complexity

As explained in Appendx A, the Cayley-Menger determinant is the determinant of an $(m+2)\times(m+2)$ symmetric matrix that relates the distances among the $m+1$ points in a set, $\Theta_{\ell}\in\mathbb{R}^{m}$, to the volume of their convex hull. Thus, the dimension of such determinants only depends on the dimension of the corresponding Euclidean space. Considering our localization algorithm in $\mathbb{R}^{2}$, when a robot finds enough (at least $m+1=3$) neighbors, it has to calculate $m+2=4$ Caley-Menger determinants of $4\times 4$ matrices to perform the inclusion test. Similarly, in $\mathbb{R}^{3}$, when a robot finds at least $m+1=4$ neighbors, it needs to calculate $m+2=5$ determinants of $5\times 5$ matrices to perform the inclusion test. Note that to perform the inclusion test in $\mathbb{R}^{2}$ and $\mathbb{R}^{3}$, a robot has to compute four areas and five volumes, respectively. Since the complexity of the computation of an $n\times n$ determinant is $\mathcal{O}(n!)$, considering at most one update per iteration, the computation complexity of the algorithm for an updating robot in $\mathbb{R}^{m}$ is $(m+2)\mathcal{O}((m+2)!)$, where $m\leq 3$.

VIII-C Communication loss and drops

In the event of a temporary communication loss/drop, a robot will not be able to perform an update and/or serve as a neighboring node in any triangulation set. However, since the algorithm is opportunistic and fully distributed, this is not going to affect the performance of the rest of the team except that it reduces the chance for the other robots to find a triangulation set. As the size of the network grows, such effect becomes negligible. However, if the failure is permanent, i.e., one robot gets isolated from the network by the end of the localizing process, it will not be able to find its location and the size of the network gets reduced by one.

VIII-D Challenges and future work

We note that a successful implementation of this algorithm requires a knowledge of the associated parameters, $\alpha,\alpha_{k},\epsilon$. We provide some insight into these questions in the simulation section of this paper and have discussed the theoretical relevance of each one throughout the paper. In future, we will study how the changes in these parameters effect the convergence or performance of the underlying localization algorithm, which may require a large number of practical experiments. As a specific example, if the statistics of the distance/motion noise variables are known, one may be able to design $\epsilon$ so the inclusion test is correct (with a very high probability) regardless of the noise parameters. As part of our future work, we are also planning to test our localization algorithm in our lab both on UGV’s (in R2) and UAV’s (in R3).

IX Conclusions

In this paper, we provide a linear distributed algorithm to localize an arbitrary number of mobile robots moving in a bounded region. We assume that each robot can measure a noisy version of its motion as well as its distance to the neighboring nodes. We consider an opportunistic algorithm such that the robots update their location estimates in ${\mathbb{R}}^{m}$, as a linear-convex combination of their $m+1$ neighbors if they lie inside their convex hull. The updating robot uses the opportunistic information exchange to refine its location estimate by using a barycentric-based convex update rule. We abstract the algorithm as an LTV system with (sub-)stochastic matrices, and show that it converges to the true robot locations under some mild regularity conditions on update weights. We also relate the dimension of motion in the network to the number of beacons required, and show that a network of mobile robots with full degrees of freedom in the motion can be localized precisely, as long as there is at least one beacon in the network. We evaluate the performance of the algorithm in the presence of noise and provide modifications to the proposed algorithm to address the undesirable effects of noise.

Appendix A Cayley-Menger determinant

Consider a set of $m+1$ points, $\Theta_{\ell}$, in $m$-dimensional Euclidean space. Cayley-Menger (CM) determinant is the determinant of an $(m+2)\times(m+2)$ symmetric matrix, which uses the pairwise distance information in the set $\Theta_{\ell}$, to compute the hypervolume, $A_{\Theta_{\ell}}$, of their convex hull, $\mathcal{C}(\Theta_{\ell})$. The CM determinant is given by

[TABLE]

in which ${\bf{1}}_{m+1}$ denotes an $m+1$-dimensional column vector of $1$s, ${\bf{D}}=\{d_{{\ell}j}^{2}\},~{}\ell,j\in\Theta_{l}$, is the $(m+1)\times(m+1)$ matrix of squared distances, $d_{{\ell}j}$, within the set, $\Theta_{\ell}$, and

[TABLE]

The second and third coefficients in the above sequence that are relevant in $\mathbb{R}^{2}$ and $\mathbb{R}^{3}$ are $-16$ and $288$, respectively.