Achieving the Desired Dynamic Behavior in Multi-Robot Systems Interacting with the Environment

TL;DR

This paper presents a methodology for controlling multi-robot systems to achieve desired dynamic behaviors during environmental interactions, ensuring safety through passivity preservation and validated via large-scale simulations.

Contribution

It introduces a local scaling approach to inter-robot coupling that guarantees passivity and desired interaction behavior in multi-robot systems.

Findings

Method guarantees passivity preservation during interactions.

Simulation results demonstrate effectiveness in large-scale systems.

Achieves desired dynamic behavior safely during environmental contact.

Abstract

In this paper we consider the problem of controlling the dynamic behavior of a multi-robot system while interacting with the environment. In particular, we propose a general methodology that, by means of locally scaling inter-robot coupling relationships, leads to achieving a desired interactive behavior. The proposed method is shown to guarantee passivity preservation, which ensures a safe interaction. The performance of the proposed methodology is evaluated in simulation, over large-scale multi-robot systems.

Equations

$$\begin{array}[l]{lllll}F^{e}=\begin{pmatrix}F^{e}_{1}\\ \vdots\\ F^{e}_{N}\end{pmatrix}\in\mathbb{R}^{3N}&F^{c}=\begin{pmatrix}F^{c}_{1}\\ \vdots\\ F^{c}_{N}\end{pmatrix}\in\mathbb{R}^{3N}&p=\begin{pmatrix}p_{1}\\ \vdots\\ p_{N}\end{pmatrix}\in\mathbb{R}^{3N}&\chi=\begin{pmatrix}\chi_{1}\\ \vdots\\ \chi_{\bar{N}}\end{pmatrix}=\begin{pmatrix}x_{1,2}\\ \vdots\\ x_{N-1,N}\end{pmatrix}\in\mathbb{R}^{3\bar{N}}&\beta=\begin{pmatrix}\beta_{1}\\ \vdots\\ \beta_{\bar{N}}\end{pmatrix}=\begin{pmatrix}\beta_{1,2}\\ \vdots\\ \beta_{N-1,N}\end{pmatrix}\in\mathbb{R}^{3\bar{N}}\end{array}$$

$$m_i{\ddot{x}}_i=w_{i}i=1,\dots ,N$$

$$w_i=-\sum _{j=1,j\ne i}^N\nabla V\left(x_{i,j}\right)-\sum _{j=1,j\ne i}^N{\beta }_{i,j}\left({\dot{x}}_i-{\dot{x}}_j\right)+F_i^c-b_i{\dot{x}}_i+F_i^e$$

$${\mathsf{N}}_i=\left\{j\ne i\text{ such that }\Vert x_{i,j}\Vert \le R\right\}$$

$$\beta_{i,j}=\left\{\begin{array}[]{ll}\beta_{i,j}\geq 0&\text{ if }\,j\in\mathsf{N}_{i}\\ 0&\text{otherwise}\end{array}\right.$$

$$m_i{\ddot{x}}_i+b_i{\dot{x}}_i+\sum _{j=1,j\ne i}^N\nabla V_{i,j}+\sum _{j=1,j\ne i}^N{\beta }_{i,j}\left({\dot{x}}_i-{\dot{x}}_j\right)=F_i^e+F_i^c$$

$$B={\mathcal{L}}_{\beta }+\mathrm{diag}\left(b_1,\dots ,b_N\right)$$

$$\left\{\begin{array}[l]{ll}\begin{pmatrix}\dot{p}\\ \dot{\chi}\end{pmatrix}=\left[\begin{pmatrix}\mathbb{O}&\mathcal{I}\\ -\cal I^{T}&\mathbb{O}\end{pmatrix}-\begin{pmatrix}B&\mathbb{O}\\ \mathbb{O}&\mathbb{O}\end{pmatrix}\right]\begin{pmatrix}{\frac{\partial H}{\partial p}}\\ {\frac{\partial H}{\partial\chi}}\end{pmatrix}+G(F^{e}+F^{c})\\ v=G^{T}\begin{pmatrix}{\frac{\partial H}{\partial p}}\\ {\frac{\partial H}{\partial\chi}}\end{pmatrix}\end{array}\right.$$

$$H=\sum _{i=1}^N{\mathcal{K}}_i(p_i)+\sum _{i=1}^N\sum _{i=j,j\ne i}^N{V}_{i,j}\ge 0$$

$$H=\sum _{i=1}^N{\mathcal{K}}_i(p_i)+\sum _{k=1}^{\bar{N}}{V}_k\ge 0$$

$$\begin{array}[l]{ll}\dot{H}=\begin{pmatrix}{\frac{\partial^{T}H}{\partial p}}&{\frac{\partial^{T}H}{\partial\chi}}\end{pmatrix}\begin{pmatrix}\dot{p}\\ \dot{\chi}\end{pmatrix}=\begin{pmatrix}{\frac{\partial^{T}H}{\partial p}}&{\frac{\partial^{T}H}{\partial\chi}}\end{pmatrix}\left[\begin{pmatrix}0&\mathcal{I}\\ -\cal I^{T}&0\end{pmatrix}-\right.\\ \left.\begin{pmatrix}B&0\\ 0&0\end{pmatrix}\right]\begin{pmatrix}{\frac{\partial H}{\partial p}}\\ {\frac{\partial H}{\partial\chi}}\end{pmatrix}+\begin{pmatrix}{\frac{\partial^{T}H}{\partial p}}&{\frac{\partial^{T}H}{\partial\chi}}\end{pmatrix}G(F^{e}+F^{c})=\\ =-{\frac{\partial^{T}H}{\partial p}}B{\frac{\partial H}{\partial p}}+(F^{e}+F^{c})^{T}v\leq(F^{e}+F^{c})^{T}v\end{array}$$

$${\int }_0^T(F^c+F^e{)}^{T}vd\tau \ge H(t)-H(0)\ge -H(0)$$

$$\alpha_{k}=\left\{\begin{array}[l]{lll}\alpha(t)>0&\mbox{if }\mathcal{I}_{\mathcal{G}}\left[i,k\right]\neq 0\text{ and }\left\|\chi_{k}\right\|\leq R\\ 1&\mbox{otherwise }\end{array}\right.$$

$$c_{k}=\left\{\begin{array}[l]{lll}\gamma(t)>0&\mbox{if }\mathcal{I}_{\mathcal{G}}\left[i,k\right]\neq 0\text{ and }\left\|\chi_{k}\right\|\leq R\\ 1&\mbox{otherwise }\end{array}\right.$$

$$\begin{gathered} ({\mathcal{I}}_{\mathcal{G}}{A}_{\mathcal{G}}(t))\otimes {\mathbb{I}}_3=({\mathcal{I}}_{\mathcal{G}}{A}_{\mathcal{G}}(t))\otimes ({\mathbb{I}}_3{\mathbb{I}}_3)= \\ =({\mathcal{I}}_{\mathcal{G}}\otimes {\mathbb{I}}_3)(A_{\mathcal{G}}(t)\otimes {\mathbb{I}}_3)=\mathcal{I}A(t) \end{gathered}$$

$$\underline{B}={\bar{\mathcal{L}}}_{\beta }+\mathrm{diag}\left(b_1,\dots ,b_N\right)$$

$${\bar{\mathcal{L}}}_{\beta }={\mathcal{I}}_{\mathcal{G}}\left(C_{\mathcal{G}}\bar{B}\right){\mathcal{I}}_{\mathcal{G}}^T$$

$$\left\{\begin{array}[l]{ll}\begin{pmatrix}\dot{p}\\ \dot{\chi}\end{pmatrix}=\left[\begin{pmatrix}\mathbb{O}&\mathcal{I}A(t)\\ -\cal I^{T}&\mathbb{O}\end{pmatrix}-\begin{pmatrix}\underline{B}&\mathbb{O}\\ \mathbb{O}&\mathbb{O}\end{pmatrix}\right]\begin{pmatrix}{\frac{\partial H}{\partial p}}\\ {\frac{\partial H}{\partial\chi}}\end{pmatrix}+G(F^{e}+F^{c})\\ v=G^{T}\begin{pmatrix}{\frac{\partial H}{\partial p}}\\ {\frac{\partial H}{\partial\chi}}\end{pmatrix}\end{array}\right.$$

$$\left\{\begin{array}[l]{ll}\begin{pmatrix}\dot{p}\\ \dot{\chi}\end{pmatrix}=\left[\begin{pmatrix}\mathbb{O}&\mathcal{I}\\ -\cal I^{T}&\mathbb{O}\end{pmatrix}-\begin{pmatrix}\underline{B}&\mathbb{O}\\ \mathbb{O}&\mathbb{O}\end{pmatrix}\right]\begin{pmatrix}{\frac{\partial H}{\partial p}}\\ A(t){\frac{\partial H}{\partial\chi}}\end{pmatrix}+G(F^{e}+F^{c})\\ v=G^{T}\begin{pmatrix}{\frac{\partial H}{\partial p}}\\ {\frac{\partial H}{\partial\chi}}\end{pmatrix}\end{array}\right.$$

$$H_s=\sum _{i=1}^N{\mathcal{K}}_i(p_i)+\sum _{k=1}^{\bar{N}}{\alpha }_k{V}_k\ge 0$$

$$H_s(t)=\sum _{i=1}^N{\mathcal{K}}_i(p_i)+\sum _{k=1}^{\bar{N}}{\alpha }_k(t)V_k$$

$${\dot{H}}_s=(F^c+F^e{)}^{T}v-\frac{{\partial }^T{H}_s}{\partial p}\underline{B}\frac{\partial H_s}{\partial p}+\frac{\partial H}{\partial t}$$

$$H_s(t)-H_s(0)\le {\int }_0^t(F^c+F^e{)}^{T}vd\tau +{\int }_0^t\frac{\partial H}{\partial \tau }d\tau$$

$$\begin{array}[l]{ll}\displaystyle\int_{0}^{t}{\frac{\partial H}{\partial\tau}}d\tau=\!\!\int_{0}^{t}\sum_{k=1}^{\bar{N}}\dot{\alpha}_{k}(\tau)V_{k}(\tau)d\tau=\\ =\!\displaystyle\sum_{k=1}^{\bar{N}}\alpha_{k}(t)V_{k}(t)-\!\sum_{k=1}^{\bar{N}}\alpha_{k}(0)V_{k}(0)-\!\int_{0}^{t}\!\sum_{k=1}^{\bar{N}}\alpha_{k}(\tau)\dot{V}_{k}(\tau)d\tau\end{array}$$

$${\Phi }_k^{+}=\left\{t\in \left[0,T\right]\text{ such that }{\dot{V}}_k(t)\ge 0\right\}$$

$${\Phi }_k^{-}=\left\{t\in \left[0,T\right]\text{ such that }{\dot{V}}_k(t)<0\right\}$$

$${\Phi }_k^{+}=\left[{\underline{t}}_{k,1}^{+},{\bar{t}}_{k,1}^{+}\right]\bigcup \dots \bigcup \left[{\underline{t}}_{k,T_k^{+}}^{+},{\bar{t}}_{k,T_k^{+}}^{+}\right]$$

$${\Phi }_k^{-}=\left[{\underline{t}}_{k,1}^{-},{\bar{t}}_{k,1}^{-}\right]\bigcup \dots \bigcup \left[{\underline{t}}_{k,T_k^{-}}^{-},{\bar{t}}_{k,T_k^{-}}^{-}\right]$$

$$0\le {\underline{t}}_{k,1}^{+}<{\bar{t}}_{k,1}^{+}<\dots <{\underline{t}}_{k,T_k^{+}}^{+}<{\bar{t}}_{k,T_k^{+}}^{+}\le T$$

$$0\le {\underline{t}}_{k,1}^{-}<{\bar{t}}_{k,1}^{-}<\dots <{\underline{t}}_{k,T_k^{-}}^{-}<{\bar{t}}_{k,T_k^{-}}^{-}\le T$$

$$\begin{array}[]{l}\displaystyle\int_{0}^{t}\sum_{k=1}^{\bar{N}}\alpha_{k}(\tau)\dot{V}_{k}(\tau)d\tau=\\ \displaystyle\sum_{k=1}^{\bar{N}}\left(\int_{\Phi_{k}^{+}}\alpha_{k}(\tau)\dot{V}_{k}(\tau)d\tau+\int_{\Phi_{k}^{-}}\alpha_{k}(\tau)\dot{V}_{k}(\tau)d\tau\right)=\\ \displaystyle\sum_{k=1}^{\bar{N}}\left(\sum_{h=1}^{T_{k}^{+}}\int_{\underline{t}_{k,h}^{+}}^{\bar{t}_{k,h}^{+}}\alpha_{k}(\tau)\dot{V}_{k}(\tau)d\tau+\sum_{h=1}^{T_{k}^{-}}\int_{\underline{t}_{k,h}^{-}}^{\bar{t}_{k,h}^{-}}\alpha_{k}(\tau)\dot{V}_{k}(\tau)d\tau\right)\end{array}$$

$$\begin{array}[]{l}-\displaystyle\int_{0}^{t}\sum_{k=1}^{\bar{N}}\alpha_{k}(\tau)\dot{V}_{k}(\tau)d\tau\leq\\ \displaystyle-\sum_{k=1}^{\bar{N}}\left(\alpha_{m}\sum_{h=1}^{T_{k}^{+}}\int_{\underline{t}_{k,h}^{+}}^{\bar{t}_{k,h}^{+}}\dot{V}_{k}(\tau)d\tau+\alpha_{M}\sum_{h=1}^{T_{k}^{-}}\int_{\underline{t}_{k,h}^{-}}^{\bar{t}_{k,h}^{-}}\dot{V}_{k}(\tau)d\tau\right)\end{array}$$

Full text

I Introduction

This paper proposes a general decentralized methodology for achieving a desired overall dynamic behavior for a multi-robot system interacting with the environment.

Typically, the behavior of a multi-robot system is defined by the interplay among basic control actions, such as aggregation, swarming, formation control, coverage and synchronization [1, 2, 3, 4]. Modifying those basic actions makes it possible to change the characteristic properties of the overall multi-robot system. Along these lines, several methods can be found in the literature that tune the inter-robot coupling actions to modify global geometric properties of the group [2, 5, 6, 7], in terms of relative positions.

Besides geometric properties, it is often of interest to regulate some topological properties of the multi-robot systems, such as connectivity [5, 8], bi-connectivity [9, 10], controllability [11] or rigidity [12].

In this paper we consider the problem of achieving a desired dynamic behavior when the multi-robot system is interacting with the environment, by means of appropriately tuning the coupling among neighboring robots. This problem was addressed in [13], where an observation and estimation scheme was defined for understanding the behavior of humans. In particular, a common scaling factor was introduced, to reduce the inter-robot forces and, thus, impose constraints on the velocities and accelerations of the robots, when needed.

However, it is worth noting that uniformly scaling down the interaction forces among all the robots might lead to a too conservative solution, where connections among robots become too loose, and the primary objective of the multi-robot system can not be correctly fulfilled.

In this paper we propose a strategy for achieving a desired interactive behavior of a multi-robot system with the environment. To this aim we will adopt a passivity based approach. In fact, guaranteeing the passivity of the multi-robot system is a sufficient condition for ensuring a stable behavior during the interaction with the, even poorly known, environment [14]. Given a general passive cooperative nominal behavior of the multi-robot system, when a robot interacts with the environment, a desired interactive behavior is achieved by nonlinearly scaling the coupling with its neighbors. The proposed control of interaction is intended as a low-level control layer, to be coupled with some nominal control action, and is designed in a local manner, in order to affect as less as possible the nominal behavior of the multi-robot system.

The paper is organized as follows. Section II introduces the notation used in the paper. Problem formulation is provided in Section III. A method for tuning the coupling gains while preserving passivity is described in Section IV. This methodology is exploited in Section V for achieving the desired viscoelastic dynamic behavior. Simulations are described in Section VI, and concluding remarks are given in Section VII.

II Notation

The symbol $\mathbb{I}_{m}\in\mathbb{R}^{m\times m}$ will indicate the identity matrix of dimension $m$, and the symbol $\mathbb{O}_{m,n}\in\mathbb{R}^{m\times n}$ will indicate the null matrix of dimension $m\times n$. For ease of notation, we will omit the dimension of the matrices when they appear clearly from the context.

Let $\Omega\in\mathbb{R}^{\rho\times\sigma}$ be a generic matrix. Then, we define $\Omega\left[i,j\right]\in\mathbb{R}$ as the element $\left(i,j\right)$ of $\Omega$. Moreover, the symbol $\otimes$ will be used to represent the Kronecker product.

III Problem Formulation

Consider a system composed of $N$ robots moving in a three-dimensional environment, whose dynamics are modeled as follows:

[TABLE]

where $x_{i}\in\mathbb{R}^{3}$ is the $i$-th robot’s position, $m_{i}>0$ is the $i$-th robot’s mass, and $w_{i}\in\mathbb{R}^{3}$ collects control inputs and all the external forces each robot is subject to.

We consider the case where each robot is controlled in such a way that some desired cooperative behavior is achieved, while robots can interact with the environment. Hence, we consider the following generic input model:

[TABLE]

where $x_{i,j}=x_{i}-x_{j}$. The terms in (2) are defined as follows.

The term $-\mkern-15.0mu\displaystyle\sum_{j=1,j\neq i}^{N}\mkern-10.0mu\nabla V\left(x_{i,j}\right)$ represents the coupling among robots. In particular, we consider each robot to interact with its neighbors, namely those robots whose distance is smaller than a certain threshold $R>0$, implementing a gradient descent of the artificial potential field $V\left(x_{i,j}\right)\geq 0$ [4, 15, 14, 16], that has a global minimum at the desired inter-robot distance $\left\|x_{i,j}\right\|=\delta_{d}>0$. The potential field is then designed in such a way that an attractive force is generated if $\delta_{d}\leq\left\|x_{i,j}\right\|\leq R$, and a repulsive force is generated if $\left\|x_{i,j}\right\|<\delta_{d}$, such that the inter-robot distance does never go below the safety value $\delta_{s}$, with $0<\delta_{s}<\delta_{d}$. A zero force is generated if two robots are too far away from each other, namely if $\left\|x_{i,j}\right\|>R$.

This kind of coupling represents an elastic interconnection among the robots, by means of nonlinear springs.

According to this definition of the coupling potential, we define the $i$-th robot’s neighborhood as follows:

[TABLE] 2. 2.

Inter-robot damping is represented by the term $-\mkern-15.0mu\displaystyle\sum_{j=1,j\neq i}^{N}\mkern-10.0mu\beta_{i,j}\left(\dot{x}_{i}-\dot{x}_{j}\right)$, where $\beta_{i,j}$ is defined as follows:

[TABLE]

Together with the definition of the previously introduced coupling term, the design of $\beta_{i,j}$ leads to defining the overall desired behavior of the multi-robot system. 3. 3.

The term $F_{i}^{c}$ represents an additional control input for the $i$-th robot, that can be utilized for achieving different objectives, such as imposing an offset [2] or obtaining complex behaviors [17]. 4. 4.

The local damping term $-b_{i}\dot{x}_{i}$, with $b_{i}>0$, represents both the viscous friction that characterizes the system and any additional damping injection obtained through a local control action. 5. 5.

The term $F_{i}^{e}$ represents the interaction force of the $i$-th robot with the environment. It can be either a real contact force, measured by means of force sensors, or a virtual force, generated by an obstacle avoidance artificial potential field [15, 18].

For ease of notation, we will hereafter define $V_{i,j}=V\left(x_{i,j}\right)$. Hence, considering the input defined in (2), the dynamics of the $i$-th robot introduced in (1) can be rewritten as follows:

[TABLE]

Besides defining how the robots coordinate among each other, (5) defines also how the multi-robot system interacts with the environment: the coupling forces among the robots define the overall viscoelastic behavior of the multi-robot system.

In this paper, we address the following problem:

Problem**.**

*Define the coupling forces among the robots in such a way that the overall multi-robot system interacts with the environment with some desired viscoelastic behavior, while preserving its overall passivity. *

IV Tuning of the Coupling Among the Robots While Preserving Passivity

In this Section we will introduce a methodology for tuning the coupling among the robots while preserving passivity. For this purpose, we will rewrite the model of the multi-robot system in port-Hamiltonian form. Define then $p_{i}=m\dot{x}_{i}$ as the $i$-th robot’s momentum, let $\bar{N}=\frac{N(N-1)}{2}$, and consider the quantities defined in (6).

Furthermore, let $\mathcal{I}_{\mathcal{G}}\in\mathbb{R}^{N\times\bar{N}}$ be the incidence matrix of the complete graph among the robots111By complete graph we refer to an undirected graph, in which each robot is represented by a node, and an edge exists among each pair of nodes.. Define also $\bar{B}=\operatorname{diag}\left(\beta\right)$. The inter-agent damping term can then be modeled utilizing the weighted Laplacian matrix $\mathcal{L}_{\beta}\in\mathbb{R}^{N\times N}$ defined, as shown in [5], as $\mathcal{L}_{\beta}=\mathcal{I}_{\mathcal{G}}\,\bar{B}\,\mathcal{I}_{\mathcal{G}}^{T}$. Define now $M=\operatorname{diag}\left(m_{1},\dots,m_{N}\right)$ and

[TABLE]

as the inertia and damping matrix of the multi-robot system, respectively. Moreover, define $\mathcal{I}=\mathcal{I}_{\mathcal{G}}\otimes\mathbb{I}_{3}$.

The model of the multi-robot system can then be given in port-Hamiltonian form as follows:

[TABLE]

where $v=\left(\dot{x}_{1}^{T}\ldots\dot{x}_{N}^{T}\right)^{T}\in\mathbb{R}^{3N}$ is the velocity vector, $G=\left(\mathbb{I}_{3N}\quad\mathbb{O}_{3\bar{N},3N}\right)^{T}$, and $H$ is the total energy of the system given by:

[TABLE]

where $\mathcal{K}_{i}(p_{i})=\nicefrac{{p_{i}^{T}p_{i}}}{{2m_{i}}}$ is the kinetic energy associated to robot $i$. For ease of notation, define $V_{k}=V\left(\chi_{k}\right)$. Hence, it is possible to rewrite (9) as:

[TABLE]

The following result can be trivially derived from [14, Proposition 1].

Proposition 1**.**

*Consider the dynamics of the multi-robot system described in port-Hamiltonian form in (8), and consider the total energy of the system given in (10). Then, the system is passive with respect to the pair $(F^{c}+F^{e},v)$. *

Proof.

Consider the definition of the total energy of the system $H$ given in (10). Then, considering the dynamics of the multi-robot system given in (8), the time derivative of $H$ can be computed as follows:

[TABLE]

Thus,

[TABLE]

which completes the proof.

Thus, the multi-robot system can safely interact with the environment, but its dynamic behavior is determined by the inter-robot coupling. We will hereafter define a methodology for scaling the coupling forces among the robots, with the objective of achieving some desired dynamic viscoelastic behavior.

For ease of notation, we will hereafter make the following assumption:

Assumption 1**.**

*For every time $t>0$, only one robot $i=1,\ldots,N$ exists such that the interaction force with the environment $F_{i}^{e}$ is different from zero. *

Without loss of generality, in the rest of the paper we will always assume the $i$-th robot to be in contact with the environment at a given time $t$.

In order to modify the dynamic viscoelastic behavior in the interaction with the environment, we need to let the $i$-th robot adjust the force that couples it with its neighbors: this is possible by introducing a scaling factor. It is worth noting that such a scaling affects the relation between the robots and the coupling actions (elastic and damping forces).

Define now $A_{\mathcal{G}}=\operatorname{diag}(\alpha_{1},\dots,\alpha_{\bar{N}})$ as a diagonal matrix containing scaling factors, one per each pair of robots. If robot $i$ is interacting with the environment, it is sufficient to scale its interaction with the neighbors. Considering the definition of the $i$-th robot’s neighborhood given in (3), the elements of $A_{\mathcal{G}}$ can be defined as follows:

[TABLE]

In a similar manner, it is possible to introduce scaling factors for modulating the inter-agent damping. Define then $C_{\mathcal{G}}=\operatorname{diag}\left(c_{1},\ldots,c_{\bar{N}}\right)$ as a diagonal matrix containing a scaling factor per each pair of robots. As in the case of $A_{\mathcal{G}}$, define the elements of $C_{\mathcal{G}}$ as follows:

[TABLE]

In this way, only the edges transmitting coupling forces between $i$ and its neighbors are scaled. In particular, the matrix transmitting the scaled elastic forces to the robots is given by $(\mathcal{I}_{\mathcal{G}}A_{\mathcal{G}}(t))\otimes\mathbb{I}_{3}$. Exploiting the properties of the Kronecker product we have that:

[TABLE]

where $A(t)=A_{\mathcal{G}}(t)\otimes\mathbb{I}_{3}$. Moreover, the matrix transmitting the scaled viscous forces to the robots is given by

[TABLE]

where

[TABLE]

Thus, the model of the multi-robot system (8) is modified as follows, with the introduction of the scaled coupling forces:

[TABLE]

Comparing (18) with (8), it is possible to note that scaling of the inter-agent damping does not change the structure of the model: in fact, both $\underline{B}$ and $B$ are positive definite: therefore, their role in the proof of Proposition 1 is analogous.

Conversely, due to the presence of the scaling matrix $A(t)$, the interconnection matrix in (18) is not skew-symmetric as in (8). As a consequence, the result of Proposition 1 can not be applied in this case. Intuitively, this is due to the fact that scaling the exchanged forces destroys the power balance among the robots. Nevertheless, we will hereafter show that passivity can still be guaranteed.

Considering time varying gains, (18) can be rewritten as

[TABLE]

As a consequence, the energy function (10) is now modified as follows:

[TABLE]

Proposition 2**.**

*Consider the dynamics of the multi-robot system described in port-Hamiltonian form in (19), and consider the total energy of the system given in (20). Then, if $\exists\alpha_{m},\alpha_{M}\in\mathbb{R}$, $0<\alpha_{m}<\alpha_{M}$, such that $\alpha_{m}\leq\alpha_{k}(t)\leq\alpha_{M}$ for any time $t\geq 0$, then the system is passive with respect to the pair $(F^{c}+F^{e},v)$. *

Proof.

Consider the following energy function:

[TABLE]

Following the same steps taken in the proof of Proposition 1, we have that

[TABLE]

and consequently, by integrating and by reminding that $\underline{B}$ is positive definite:

[TABLE]

Integrating by parts we get

[TABLE]

For any time interval $\left[0,T\right]$, and for any $k=1,\ldots,\bar{N}$, it is possible to define the following subsets:

[TABLE]

The subsets $\Phi_{k}^{+}$ and $\Phi_{k}^{-}$ are the union of a finite number of disjoint time intervals222Pathological situations might exist in which the number of time intervals that compose $\Phi_{k}^{+}$ and $\Phi_{k}^{-}$ is not finite. However, in practical situations this does never happen.. Let $T_{k}^{+}$ and $T_{k}^{-}$ be the number of time intervals that compose $\Phi_{k}^{+}$ and $\Phi_{k}^{-}$, respectively. Let $\underline{t}_{k,h}^{+},\bar{t}_{k,h}^{+}$ be the initial and final times of the $h$-th time interval of $\Phi_{k}^{+}$, and let $\underline{t}_{k,h}^{-},\bar{t}_{k,h}^{-}$ be the initial and final times of the $h$-th time interval of $\Phi_{k}^{-}$. Then, the subsets can be defined as follows:

[TABLE]

where

[TABLE]

Hence, we can rewrite the integral in (24) as follows:

[TABLE]

Since $0\leq\alpha_{m}\leq\alpha(t)\leq\alpha_{M}$, according to the definition given in (26), it is possible to obtain the following inequality:

[TABLE]

Hence, from (24) we obtain the following:

[TABLE]

Thus, from (23), the following inequality can be derived:

[TABLE]

Hence, considering the definition of $H_{s}(t)$ given in (21), it is possible to obtain the following inequality:

[TABLE]

Since both $\alpha_{m}$ and $\alpha_{M}$ are positive, and since both the kinetic energy $K_{i}\left(\cdot\right)$ and the potentials $V_{k}\left(\cdot\right)$ are positive, it is possible to obtain the following:

[TABLE]

which proves the passivity.

V Tunable Interaction with the Environment

In this Section we will show how to utilize the methodology introduced so far for adjusting the parameters of the inter-robot coupling, in a local manner, in order to achieve a desired viscoelastic dynamic behavior in the interaction with the environment. We will hereafter assume that the robots’ mass is sufficiently small, such that inertial forces can be neglected. Hence, the force robot $i$ applies to the environment is equal to

[TABLE]

Since the elastic coupling term between any two robots is only a function of their relative positions, it is always possible to write it as follows:

[TABLE]

For ease of notation, we will hereafter omit the dependency of $\kappa_{i,j}$ on $x_{i},x_{j}$. Hence, (34) can be rewritten as follows:

[TABLE]

This force can be modeled as a single standard viscoelastic force as follows:

[TABLE]

Let $\Delta$ be a constant representing the rest-length of a standard elastic element, defined based on the application. Then, the nominal stiffness $\kappa_{n}$ and rest position $\bar{x}$ of the spring are defined as follows:

[TABLE]

[TABLE]

The value of $\bar{x}$ is then defined according to (38) in such a way that the following holds:

[TABLE]

Parameters $\bar{v}$ and $\beta_{n}$ are defined in an analogous manner.

It is worth noting that $\bar{x}$ plays the role of the desired position for robot $i$ in stiffness control. Hence, it is possible to define a desired dynamic behavior for the multi-robot system, in terms of a desired viscoelastic dynamics, as follows:

[TABLE]

for some desired $\kappa_{d}>0$, $\beta_{d}>0$. As discussed in Section IV, the damping coefficient can be freely adjusted with an appropriate choice of $\gamma\left(t\right)$. Therefore, imposing $\beta_{n}=\beta_{d}$, we obtain

[TABLE]

Conversely, the desired elastic stiffness can be achieved exploiting the results of Proposition 2.

In particular, it is possible to tune the coupling between the $i$-th robot and its neighbors utilizing the parameter $\alpha(t)$: the objective is that of minimizing the difference between $F_{i}$ and $F_{d}$. Namely, considering the elastic terms in (37) and (41), it is necessary to minimize the following cost function:

[TABLE]

Hence, we can define the following simple quadratic optimization problem:

[TABLE]

where $\alpha_{M}>\alpha_{m}>0$ are the upper- and lower-bounds for $\alpha(t)$, respectively, defined according to Proposition 2.

In a decentralized multi-robot system, decisions are taken in a local manner, without any centralized elaboration unit. As a consequence, if the $i$-th robot is in contact with the environment, it will locally solve the optimization problem in (44) and find the desired value $\alpha^{\star}$. It is then necessary for the $i$-th robot to broadcast this value to its neighbors, in such a way that the coupling actions can be tuned as required.

This procedure needs to be performed as soon as the $i$-th robot perceives an interaction force with the environment. Furthermore, the desired value $\alpha^{\star}$ needs to be periodically recomputed, every $\bar{T}>0$ seconds, based on the current force measurements.

VI Simulations

In this Section we describe the results of the evaluation of the proposed control method. In particular, several simulations were performed in an environment developed in MATLAB®. Specifically, a variable number of three-dimensional double integrator robots, modeled according to (1), was considered. Let the three-dimensional environment be defined by the $\left(\mathsf{x,y,z}\right)$ axes.

We utilized the coupling artificial potential field defined in [15], with the following parameter set: $\delta_{s}=5$, $\delta_{d}=15$, $R=22$. The additional control input $F_{i}^{c}\in\mathbb{R}^{3}$ was utilized for imposing a motion of the group of robots along the $\mathsf{x}$ axis. A point obstacle was then placed in the environment, and a repulsive artificial potential field was activated for those robots whose distance from the obstacle was smaller than $\delta_{d}$, in order to guarantee that the distance remained larger than $\delta_{s}$.

For each simulation run, the number $N$ of robots was defined, and initial positions were randomly chosen. Since, as detailed in Section V, the parameter $\gamma\left(t\right)$ can be arbitrarily chosen to achieve the desired damping, the evaluation focused on the stiffness tuning. Then, two cases were considered:

The nominal case, where we utilized a constant coupling gain $\alpha\left(t\right)=30$. 2. 2.

The tunable stiffness case, where we utilized the desired stiffness $\kappa_{d}=1$. In this case, the optimization problem (44) was defined with $\alpha_{m}=10^{-4}$, and $\alpha_{M}=10^{2}$.

The results of a representative simulation run are depicted in Fig. 1, where we utilized, $N=64$ robots: the tunable stiffness case is compared to the nominal one.

The value of the cost function $f\left(\alpha\right)$ defined in (43) is depicted in Fig. 1(a). The cost function was evaluated only when a robot was in contact with the obstacle. As expected, the value of $f\left(\alpha\right)$ is typically very large for the nominal case, while it becomes very small for the tunable stiffness case.

In order to evaluate the deviation from the nominal behavior, we measured the position of the barycenter of the multi-robot system. The percentage deviation in the tunable stiffness case with respect to the nominal one is depicted in Fig. 1(b) where, due to space limitations, only the component along the $\mathsf{x}$ axis is depicted. It is possible to note that large deviations are observed only during limited periods of time, in particular when a robot is in contact with the obstacle.

The accompanying video shows a few examples of simulation runs, where different numbers of robots were utilized. It is possible to note that, in the nominal case, the multi-robot system behaves as a rigid body, while in the tunable stiffness case the coupling between the robots is locally weakened, when in contact with the obstacle, in order to achieve the desired dynamic behavior.

VII Conclusions

In this paper we addressed the problem of controlling the dynamic behavior of a multi-robot system interacting with the environment. In particular, we proposed a general methodology that, introducing a local scale factor on the inter-robot couplings, leads to achieving a desired overall viscoelastic dynamic behavior. The proposed method was shown to guarantee passivity preservation, thus ensuring a safe interaction.

Throughout the paper, we assumed only one robot at a time to be in contact with the environment. Future work will aim at relaxing this assumption, considering multiple robots simultaneously in contact with the environment.