Formation Control and Obstacle Avoidance for Multiple Robots Subject to Wheel-Slip

Abstract

An adaptive formation control law for non-holonomic dynamic robots based on an artificial potential function method in the presence of lateral slip and parametric uncertainties is presented to organize multiple robots into formation. It is formulated to achieve the smooth control of the translational and rotational motion of a group of mobile robots while keeping a prescribed formation and avoiding inter-robot and obstacle collisions. In order to improve the formation control method effectively and reduce the distortion shape, the virtual leader-following method is proposed for every robot and an improved optimal assignment algorithm is used to solve multi-targets optimal assignment for the formation problem. Simulation results are provided to validate the theoretical results.

1. Introduction

Multi-robot systems (MRSs) have important advantages over single robot systems, such as flexibility, robustness, efficiency and redundancy. In several applications, a group of robots are required to accomplish the desired tasks cooperatively, maintaining a specific pattern called a ‘formation’. The formation control (FC) of a MRS is a technology for keeping the robots in the scheduled formation and adapting to the environment while moving to the target position [1]. The formation control of multiple robots is very important in a number of areas, such as: (a) personal robotics systems where small independent robots can carry out the task of moving an object, such as a table, in a cooperative manner in those situations were a single robot cannot perform the task alone; (b) military rescue and reconnaissance operations, where vehicles are often required to maintain a close formation; (c) the self-reconfiguring of metamorphic robotic systems; (d) controlling the formation of satellites for scientific and planetary exploration; and (e) autonomous vehicle highway systems where it is necessary to control a platoon of cars in formation while allowing them to perform such manoeuvres as lane changes, merger, or obstacle avoidance. In many other applications, the use of MRS rather than a single robot is urgent in improving redundancy, sensitivity and robustness as well as dealing with dangerous or distributed tasks. The study of formation by MRS extends the classical problems of single robots through new issues, like motion planning, task decomposition and multi-task assignment, searching and mapping, network communications, localization and path planning, etc[2].

Formation control is an important research field for MRS [3]. Formation control is viewed as a kind of information consensus in which robots or agents interact with each other using various sensors and communication techniques. The main goal is to coordinate a swarm of mobile agents or robots so as to achieve a desired formation pattern, avoiding inter-agent collisions at the same time.

The most fundamental aspect of formation control is the method of control used in the multi-robot or swarm system. The control method follows three different types, including (a) centralized, (b) completely decentralized and (c) hybrid. Robots need to maintain a specific formation shape by maintaining the correct formation position among other robots. Four main formation control methods have been adapted to address the formation problem of multiple agents, including the behaviour-based [4], leader-following [5], graph-theoretic [6] and virtual structure approaches [7].

In the behaviour-based [4] approach to formation control, each agent has several desired behaviours and the control action for each swarm member is defined by a weighting of the relative importance of each behaviour. In addition there are many formation control strategies which utilize artificial potential fields (APFs). Behaviour-based methods and APF are often combined in the formation control of mobile robot systems. In these approaches, each robot has basic motor schemas which generate a vector representing the desired behavioural response to sensory input. These motor schemas include behaviours such as obstacle avoidance, formation maintenance, target approaching and random wandering. However, the disadvantage of this approach is its complexity, its limited ability to maintain precise geometric formations and its vulnerability to a local minimum trap.

In most studies on formation control, only kinematic models for the non-holonomic robots have been considered [8]. The dynamics of wheeled mobile robots and the possible uncertainty within the dynamics are not considered. In many practical applications, most non-holonomic systems are dynamic systems. When one considers systems with non-holonomic kinematic and dynamic constraints, the formation control problem becomes more challenging [9,10]. A decentralized control scheme which achieves dynamic formation control and collision avoidance for a group of kinematic non-holonomic robots is represented [11]. The feedback control of a group of non-holonomic dynamic systems with uncertainty is considered based on a consensus strategy and back-stepping techniques [12]. A formation control design based on the potential function approach to the non-holonomic case is extended [8] so that the robot motion tends to follow the vector field defined by the potential function gradient while obeying the non-holonomic constraints. Pereira, A solution for the formation control of a group of non-holonomic uncertain agents and a trajectory track is proposed by this system [16]. A set of artificial potential field functions to navigate a terrain containing some obstacles for a team of non-holonomic mobile robots is proposed [14] while maintaining a desired formation and using a leader-following strategy. A formation control scheme for a group of mobile robots based on multi-objective potential force so as to reduce the global orientation error asymptotically to zero while maintaining proper formation for a target configuration is presented [15] using a translational force input and a rotational torque input. However, these papers do not consider kinematic-dynamic constraints for non-holonomic mobile robots subject to wheel-slip in the formation control problem. In [13], a complete dynamic model of a non-holonomic mobile robot in a multi-robot formation is proposed and applied to an inverse dynamics technique with uncertainties and lateral slip so as to design a formation control. In [17], it presented an approach to autonomous pursuit evasion by a wheeled mobile robot in the presence of wheeled slip using a feedback linearization controller to capture wheel-slip affection. In [20], an effect of slip on the formation control problem of wheeled mobile robots is investigated and the slip-traction characteristic is explicitly modelled and integrated into the WMR dynamics. However, this literature cannot explicitly solve robot path-planning and obstacle avoidance in the formation control design subject to wheel-slip. Among these references, [20] is the closest to our paper, but there are certain differences from the approach proposed. In [20], the main contribution is to investigate how a stable formation control algorithm with wheel-slip is introduced and thereby design a slip-based controller that can stabilize the formation of the WMR in the presence of wheel-slip. The contribution of our paper looks to investigate the impact of slip on formation control when slip is introduced into the overall dynamics of the MRS and design an adaptive formation control law for non-holonomic dynamic robots based on the artificial potential function method in the presence of lateral slip and obstacles. We use this control law to achieve the smooth control of the translational and rotational motion of a group of mobile robots while maintaining the prescribed formation and avoiding inter-robot and obstacle collisions. Another contribution of this paper is to propose a virtual target for every robot where the robot has to make an approach based on a prescribed formation geometry shape, and then utilize an optimal assignment algorithm to solve multi-target optimal assignment for the formation problem in order to avoid robots frequently changing their positions roles and so obtain a steady formation shape for MRS.

The rest of this paper is organized as follows. In Section 2, the kinematics and dynamics equations of a wheeled robot in the presence of lateral slip and dynamic uncertainties are described. Section 3 describes an artificial potential field function, including an interactive force between robots, an attraction force to target and a repulsion force to obstacles. In Section 4, a regulating control law stabilizing a team of non-holonomic robots in an unknown environment is proposed and its stability is proved. Section 5 introduces an improved Hungarian algorithm to solve multi-target optimal assignment for the formation problem. In Section 6, the evaluation method for formation shape control is discussed. Some simulation results are given in Section 7. Finally, in Section 8, the study is discussed and concluded.

2. Dynamic modelling of MRS

2.1. Formulation of the non-holonomic mobile robot system

Consider a team of n robots with double integrator dynamics. As indicated in Fig. 1(a), the number of the team is a non-holonomic wheeled mobile robot. The X-O-Y is the global coordinate system and X _C -G-Y _C is the local coordinate frame fixed on the mobile robot with the origin at the robot's centre of gravity, X _C in the heading direction of the robot and Y _C in the lateral direction. The points G, B, R, L, E, P and C represent the centre of mass of the robot, the wheel baseline centre, the ground contact points of the right and left wheels, the ground contact point of the castor, the tool or sensor location and the centre point of geometry of the robot that is required to follow a path. Forces and moments are exerted at points R, L, E and C. F _rrx and F_rry are the longitudinal and lateral tire traction forces of the right wheel, F_rlx, and F_rly are the longitudinal and lateral tire traction forces of the left wheel, F_cx and F_cy are the longitudinal and lateral tire forces exerted on C by the castor,

$\begin{array}{l}M(q)\ddot{q}+C(q,\dot{q})\dot{q}+F(\dot{q})+G(q)\\ \begin{array}{ccc}& & \end{array}+{\tau }_d=B(q)\tau -A^T(q)\lambda \end{array}$

Figure 1.

A non-holonomic mobile robot

and q ∈ ℛ^n×1 are the longitudinal and lateral tire forces exerted on E, which is the tool location, and q̇ ∈ ℛ ^n×1 is the moment exerted by the tool. M(q) ∈ ℛ^n×n and C(q,q̇) ∈ ℛ^n×n are the longitudinal and lateral velocities of the centre of the whole robot mass and B(q) ∈ ℛ^n×(n-m) and τ ∈ ℛ^(n-m)×1 are the angular velocity and heading of the robot.

In this work, we want to investigate the formation control problem of a non-holonomic mobile robot when the friction is constrained by terrain characteristics and slip occurs due to insufficient friction. As a result, we want to include the slip in the dynamics of the system. We introduce the new slip coordinate vector δ ∈ R^m×1. ξ _l and ξ _r are longitudinal slip displacements for the left and right wheel, respectively while e˜ _l and e˜ _r are lateral slip displacements. The kinematic relationships for a non-holonomic mobile robot can be described by the following equation: h_p = [x y]^T is the point that is required to track a trajectory.

The position of the robot in an inertial Cartesian frame {O,X,Y} is completely specified by the general position vector: $q={\left[\begin{array}{ccccc}{x}_P& y_P& \theta & u& w\end{array}\right]}^T$ (1)

where x_P, y_P are the coordinates of the centre of mass P_P in the world and θ is the heading angle of the robot, u, w are the longitudinal velocities of the centre of mass and the angular velocities of the robot, respectively.

$\begin{array}{l}\dot{u}=(F_{rl{x}^{\prime }}+F_{rr{x}^{\prime }}+F_{e{x}^{\prime }}+F_{c{x}^{\prime }})/m+vw\\ \dot{v}=(F_{rl{y}^{\prime }}+F_{rr{y}^{\prime }}+F_{e{y}^{\prime }}+F_{c{y}^{\prime }})/m-uw\\ \dot{w}=\frac{1}{I_z}[d(F_{rr{x}^{\prime }}-F_{rl{x}^{\prime }})-b(F_{rl{y}^{\prime }}+F_{rr{y}^{\prime }})\\ \begin{array}{cc}& \end{array}+(e-b)F_{e{y}^{\prime }}+(c-b)F_{c{y}^{\prime }}+{\tau }_e]\\ {\dot{\omega }}_l=\frac{1}{I_{lmw}}({\tau }_l-F_{rl{x}^{\prime }}{R}_{\ln t}-B_{lmw}{\omega }_l)\\ {\dot{\omega }}_r=\frac{1}{I_{rmw}}({\tau }_r-F_{rr{x}^{\prime }}{R}_{rnt}-B_{rmw}{\omega }_r)\end{array}$ (2)

where I_rmw and I_lmw are the moment of inertia of the right and left wheel system, which includes the motor, the gearbox, the encoder and the wheel, respectively. B_rmw and B_lmw are the viscous friction coefficients of the wheel system, which includes the motor, the gearbox, the encoder and wheel, respectively. R_rnt and R_lnt are the nominal radii of the tires for the robot.

The kinematics of point P is: $\begin{array}{l}{\dot{x}}_P=u\cos \theta -v\sin \theta -(a-b)w\sin \theta \\ {\dot{y}}_P=u\sin \theta +v\cos \theta +(a-b)w\cos \theta \end{array}$

where τ _l and τ _r are the torque vectors generated by the left and right motors.

$\{\begin{array}{c}{\tau }_r=k_{ar}(V_r-k_{br}{w}_r)/R_{ar}\\ {\tau }_l=k_{al}(V_l-k_{bl}{w}_l)/R_{al}\end{array}$ (3)

where V_r and V_l are the input voltages for the right and left motors, k_br and k_bl are equal to the voltage constant multiplied by the right gear ratio and the left gear ratio, respectively, R_ar is the electric resistance constant for the right motor, R_al is the electric resistance constant for the left motor and k_ar and k_al are the torque constant multiplied by the gear ratio for the right and left motors, respectively.

In [17], u, v are the longitudinal and lateral velocities of the centre of mass and w is the angular velocity of the robot: $\begin{array}{l}u=0.5\left[r(w_r+w_l)+({\dot{\xi }}_r+{\dot{\xi }}_l)\right]\\ w=\frac{1}{2d}\left[r(w_r-w_l)+({\dot{\xi }}_r-{\dot{\xi }}_l)\right]\\ v=\frac{b}{2d}\left[r(w_r-w_l)+({\dot{\xi }}_r-{\dot{\xi }}_l)\right]+0.5({\dot{\eta }}_r+{\dot{\eta }}_l)\end{array}$ (4)

The lateral and longitudinal traction forces are functions of the slip angle and slip ratio [20], respectively, and are modelled following the magic formula tyre model [21] in the literature. The traction force between a wheel and a surface is modelled as in [21]: $\begin{array}{l}{F}_{ey}(s\alpha )=D_y\sin [C_y{\tan }^{-1}\{(1-E_y)B_y(s\alpha +S_H)\\ \begin{array}{ccc}& & \end{array}+E_y{\tan }^{-1}(B_y(s\alpha +S_H))\}]+S_V\\ F_{ex}(sr)=D_x\sin [C_x{\tan }^{-1}\{(1-E_x)B_x(sr+S_H)\\ \begin{array}{ccc}& & \end{array}+E_x{\tan }^{-1}(B_x(sr+S_H))\}]+S_V\end{array}$

where B_x, By are longitudinal and lateral stiffness factors, respectively, C_x,C_y are shape factors, D_x,D_y are peak values, E_x, E_y are the curvature factor, S_H is the horizontal shift and S_V is the vertical shift. The slip ratio sr and the slip angle sα are defined as: $s{r}_i=\frac{R_{\mathrm{int}}{\omega }_i-{\dot{\xi }}_i}{{\dot{\xi }}_i},s{\alpha }_i={\tan }^{-1}(\frac{{\dot{\eta }}_i}{{\dot{\xi }}_i}),i=l,r$

In general, most commercial mobile robots have low level PID controllers to track input reference velocities (for example, for longitudinal velocity and angular velocity) instead of the motor voltage being driven directly.

Therefore, the longitudinal u_ref and angular reference velocities w_ref as control signals provide an appropriate method for robot control. In this paper, PD velocity controllers have been considered, as described by the following equation: $\left[\begin{array}{c}{V}_u\\ V_w\end{array}\right]=\left[\begin{array}{c}{k}_{PT}(u_{ref}-u_{idea})-k_{DT}{\dot{u}}_{idea}\\ k_{PR}(w_{ref}-w_{idea})-k_{DR}{\dot{w}}_{idea}\end{array}\right]$ (5)

where: $\begin{array}{l}{u}_{idea}=0.5[r(w_r+w_l)]\begin{array}{c}\end{array},w_{idea}=\frac{1}{2d}\left[r(w_r-w_l)\right]\\ V_u=0.5(V_l+V_r)\begin{array}{c}\end{array},V_w=0.5(V_r-V_l)\end{array}$ (6)

Re-arranging the third equation of (2) and the equations of (3), we obtain: $\begin{array}{l}{F}_{rl{x}^{\prime }}=\frac{k_{al}}{R_{al}{R}_{lnt}}(V_l-k_{bl}{w}_l)-\frac{B_{lmw}}{R_{lnt}}{w}_l-\frac{I_{lmw}}{R_{lnt}}{\dot{w}}_l\\ F_{rr{x}^{\prime }}=\frac{k_{ar}}{R_{ar}{R}_{rnt}}(V_r-k_{br}{w}_r)-\frac{B_{rmw}}{R_{rnt}}{w}_r-\frac{I_{rmw}}{R_{rnt}}{\dot{w}}_r\end{array}$ (7)

From equations (4), we can obtain: $\begin{array}{l}{w}_r+w_l=\frac{2u}{r}-\frac{1}{r}({\dot{\xi }}_r+{\dot{\xi }}_l)\\ w_r-w_l=\frac{2dw}{r}-\frac{1}{r}({\dot{\xi }}_r-{\dot{\xi }}_l)\end{array}$ (8)

Combining equations (8) and (6), we obtain: $\begin{array}{l}{V}_r=k_{PT}(u_{ref}-u+0.5({\dot{\xi }}_r+{\dot{\xi }}_l))-k_{DT}(\dot{u}-0.5({\ddot{\xi }}_r+{\ddot{\xi }}_l))\\ \begin{array}{cc}& \end{array}+k_{PR}(w_{ref}-w+\frac{1}{2d}({\dot{\xi }}_r-{\dot{\xi }}_l))-k_{DR}(\dot{w}-\frac{1}{d}({\ddot{\xi }}_r-{\ddot{\xi }}_l))\\ V_l=k_{PT}(u_{ref}-u+0.5({\dot{\xi }}_r+{\dot{\xi }}_l))-k_{DT}(\dot{u}-0.5({\ddot{\xi }}_r+{\ddot{\xi }}_l))\\ \begin{array}{cc}& \end{array}-k_{PR}(w_{ref}-w+\frac{1}{2d}({\dot{\xi }}_r-{\dot{\xi }}_l))+k_{DR}(\dot{w}-\frac{1}{d}({\ddot{\xi }}_r-{\ddot{\xi }}_l))\end{array}$ (9)

According to equations (3), (7) and (9), we can obtain: $\begin{array}{l}{B}_{11}\dot{w}=B_{12}{u}_{ref}+B_{13}u+B_{14}{w}_{ref}+B_{15}w+B_{16}{\dot{\xi }}_r\\ \begin{array}{cc}& \end{array}+B_{17}{\dot{\xi }}_l+B_{18}\dot{u}+B_{19}{\ddot{\xi }}_r+B_{20}{\ddot{\xi }}_l+e{F}_{e{y}^{\prime }}+c{F}_{c{y}^{\prime }}\\ \begin{array}{cc}& \end{array}+{\tau }_e+mbuw-\frac{1}{2}mb{\ddot{\eta }}_r-\frac{1}{2}mb{\ddot{\eta }}_l\end{array}$ (10)

These parameters for the above equations (10) are as follows: $\begin{array}{c}{B}_{11}=I_z-m{b}^2+\frac{d{I}_{rmw}}{r{R}_{rnt}}-\frac{d{I}_{lmw}}{r{R}_{lnt}}+\frac{k_{ar}{k}_{DR}}{R_{ar}{R}_{rnt}}-\frac{k_{al}{k}_{DR}}{R_{al}{R}_{lnt}}\\ B_{12}=\frac{d\cdot k_{PT}}{2}(\frac{k_{ar}}{R_{ar}{R}_{rnt}}-\frac{k_{al}}{R_{al}{R}_{lnt}})\\ B_{13}=\frac{d}{2}(\frac{B_{lmw}}{r{R}_{lnt}}+\frac{k_{al}{k}_{bl}}{r{R}_{al}{R}_{lnt}}-\frac{k_{ar}{k}_{br}}{r{R}_{ar}{R}_{rnt}}-\frac{B_{rwm}}{R_{rnt}})\\ \begin{array}{cc}& \end{array}-\frac{d{k}_{PT}}{2}(\frac{k_{ar}}{R_{ar}{R}_{rnt}}-\frac{k_{al}}{R_{al}{R}_{lnt}})\\ B_{14}=\frac{d\cdot k_{PR}}{2}(\frac{k_{ar}}{R_{ar}{R}_{rnt}}+\frac{k_{al}}{R_{al}{R}_{lnt}})\\ B_{15}=-\frac{d}{2}(\frac{k_{ar}{k}_{PR}}{R_{ar}{R}_{rnt}}+\frac{k_{al}{k}_{PR}}{R_{al}{R}_{lnt}})\\ +\frac{d}{2}(\frac{B_{lmw}}{r{R}_{lnt}}+\frac{k_{al}{k}_{bl}}{r{R}_{al}{R}_{lnt}}-\frac{B_{rmw}}{r{R}_{rnt}}-\frac{k_{ar}{k}_{br}}{r{R}_{ar}{R}_{rnt}})\\ B_{16}=\frac{d}{2}[\frac{k_{ar}}{2R_{ar}{R}_{rnt}}(k_{PT}+\frac{k_{PR}}{d})-\frac{k_{al}}{2R_{al}{R}_{lnt}}(k_{PT}-\frac{k_{PR}}{d})]\\ \begin{array}{cc}& \end{array}+\frac{d}{2r}[\frac{k_{ar}{k}_{br}}{R_{ar}{R}_{rnt}}+\frac{B_{rwm}}{R_{rnt}}]\\ B_{17}=\frac{d}{2}[(\frac{k_{ar}}{2R_{ar}{R}_{rnt}}-\frac{k_{al}}{2R_{al}{R}_{lnt}})(k_{PT}-\frac{k_{PR}}{d})]\\ \begin{array}{cc}& \end{array}+\frac{d}{2r}(\frac{k_{al}{k}_{bl}}{R_{al}{R}_{lnt}}+\frac{B_{lwm}}{R_{lnt}})\\ B_{18}=-\frac{d}{2}\left[\frac{k_{ar}{k}_{DT}}{R_{ar}{R}_{rnt}}-\frac{k_{al}{k}_{DT}}{R_{al}{R}_{lnt}}+\frac{I_{rmw}}{r\cdot R_{rnt}}-\frac{I_{lmw}}{r\cdot R_{\ln t}}\right]\\ B_{19}=\frac{d}{2}\left[\frac{k_{ar}}{2R_{ar}{R}_{rnt}}(k_{DT}+\frac{k_{DR}}{d})-\frac{k_{al}}{2R_{al}{R}_{lnt}}(k_{DT}-\frac{k_{DR}}{d})\right]+\frac{d\cdot I_{rmw}}{2r\cdot R_{rnt}}\\ B_{20}=\frac{d}{2}\left[\frac{k_{ar}}{2R_{ar}{R}_{rnt}}(k_{DT}-\frac{k_{DR}}{d})-\frac{k_{al}}{2R_{al}{R}_{lnt}}(k_{DT}+\frac{k_{DR}}{d})\right]-\frac{d\cdot I_{lmw}}{2r\cdot R_{lnt}}\end{array}$

From the first expression of equations (2), (7) and (9), we can obtain: $\begin{array}{l}(m-D_{15})\dot{u}=D_{11}{u}_{ref}+D_{12}u+D_{13}{\dot{\xi }}_r+D_{14}{\dot{\xi }}_l+D_{16}{\ddot{\xi }}_r+D_{17}{\ddot{\xi }}_l\\ \begin{array}{c}\end{array}+D_{18}{w}_{ref}+D_{19}w+(F_{e{x}^{\prime }}+F_{c{x}^{\prime }})+mb{w}^2+\frac{1}{2}m({\dot{\eta }}_r+{\dot{\eta }}_l)\end{array}$

Among the above equations (11), these parameters are: $\begin{array}{c}{D}_{11}=\frac{k_{ar}{k}_{PT}}{R_{ar}{R}_{rnt}}+\frac{k_{al}{k}_{PT}}{R_{al}{R}_{lnt}}\\ D_{12}=-\frac{k_{ar}{k}_{PT}}{R_{ar}{R}_{rnt}}-\frac{k_{al}{k}_{PT}}{R_{al}{R}_{lnt}}-\frac{k_{ar}{k}_{br}}{r{R}_{ar}{R}_{rnt}}\\ \begin{array}{cc}& \end{array}-\frac{B_{lmw}}{r\cdot R_{lnt}}-\frac{k_{al}{k}_{bl}}{r\cdot R_{al}{R}_{lnt}}-\frac{B_{rmw}}{r\cdot R_{rnt}}\\ D_{13}=\frac{k_{ar}{k}_{PT}}{2R_{ar}{R}_{rnt}}+\frac{k_{al}{k}_{PT}}{2R_{al}{R}_{lnt}}+\frac{k_{ar}{k}_{PR}}{2d\cdot R_{ar}{R}_{rnt}}\\ -\frac{k_{al}{k}_{DR}}{2d\cdot R_{al}{R}_{lnt}}-\frac{k_{ar}{k}_{br}}{r\cdot R_{ar}{R}_{rnt}}+\frac{B_{rmw}}{r\cdot R_{rnt}}\\ D_{14}=\frac{k_{ar}{k}_{PT}}{2R_{ar}{R}_{rnt}}+\frac{k_{al}{k}_{PT}}{2R_{al}{R}_{lnt}}-\frac{k_{ar}{k}_{PR}}{2d\cdot R_{ar}{R}_{rnt}}\\ +\frac{k_{al}{k}_{DR}}{2d\cdot R_{al}{R}_{lnt}}+\frac{B_{lmw}}{r\cdot R_{lnt}}-\frac{k_{al}{k}_{bl}}{r\cdot R_{al}{R}_{lnt}}\\ D_{15}=-\frac{k_{ar}{k}_{DT}}{R_{ar}{R}_{rnt}}-\frac{k_{al}{k}_{DT}}{R_{al}{R}_{lnt}}-\frac{I_{rmw}}{r\cdot R_{rnt}}-\frac{I_{lmw}}{r\cdot R_{lnt}}\\ D_{16}=\frac{k_{ar}{k}_{PT}}{2R_{ar}{R}_{rnt}}+\frac{k_{al}{k}_{PT}}{2R_{al}{R}_{lnt}}+\frac{I_{rmw}}{r\cdot R_{rnt}}\\ D_{17}=\frac{k_{ar}{k}_{PT}}{2R_{ar}{R}_{rnt}}+\frac{k_{al}{k}_{PT}}{2R_{al}{R}_{lnt}}+\frac{I_{lmw}}{r\cdot R_{lnt}}\\ D_{18}=\frac{k_{ar}{k}_{PR}}{R_{ar}{R}_{rnt}}-\frac{k_{al}{k}_{PR}}{R_{al}{R}_{lnt}}\\ D_{19}=-\frac{k_{ar}{k}_{PR}}{R_{ar}{R}_{rnt}}+\frac{k_{al}{k}_{PR}}{R_{al}{R}_{lnt}}-\frac{d\cdot k_{al}{k}_{bl}}{r\cdot R_{al}{R}_{lnt}}\\ \begin{array}{cc}& \end{array}+\frac{d\cdot k_{ar}{k}_{br}}{r\cdot R_{ar}{R}_{rnt}}+\frac{d\cdot B_{lmw}}{r\cdot R_{\ln t}}-\frac{d\cdot B_{rmw}}{r\cdot R_{rnt}}\end{array}$

In order to simplify dynamic modelling for a robot with slip, we assume that the following conditions hold: $k_{ar}=k_{al},R_{ar}=R_{al},R_{rnt}=R_{lnt},I_{rmw}=I_{lmw},k_{br}=k_{bl}$ (12)

The parameters for the above equations (10) are changed to: $\begin{array}{c}{B^{\prime }}_{11}=I_z-m{b}^2,{B^{\prime }}_{12}=0,{B^{\prime }}_{13}=0,{B^{\prime }}_{14}=\frac{d\cdot k_{PR}{k}_{ar}}{R_{ar}{R}_{rnt}}\\ {B^{\prime }}_{15}=-\frac{d\cdot k_{ar}{k}_{PR}}{R_{ar}{R}_{rnt}},{B^{\prime }}_{16}=\frac{1}{2}[\frac{k_{ar}{k}_{PR}}{R_{ar}{R}_{rnt}}]+\frac{d}{2r}[\frac{k_{ar}{k}_{br}}{R_{ar}{R}_{rnt}}+\frac{B_{rwm}}{R_{rnt}}]\\ {B^{\prime }}_{17}=-\frac{d}{2r}(\frac{k_{al}{k}_{bl}}{R_{al}{R}_{lnt}}+\frac{B_{lwm}}{R_{lnt}}),{B^{\prime }}_{18}=0\\ {B^{\prime }}_{19}=\frac{d\cdot I_{rmw}}{2r\cdot R_{rnt}}+\frac{1}{2}\frac{k_{ar}{k}_{DR}}{R_{al}{R}_{lnt}}\\ {B^{\prime }}_{20}=-\frac{k_{al}{k}_{DR}}{2R_{al}{R}_{lnt}}-\frac{d\cdot I_{lmw}}{2r\cdot R_{lnt}},{D^{\prime }}_{11}=\frac{2k_{ar}{k}_{PT}}{R_{ar}{R}_{rnt}}\\ {D^{\prime }}_{12}=-\frac{2k_{ar}{k}_{PT}}{R_{ar}{R}_{rnt}}-\frac{2k_{ar}{k}_{br}}{r{R}_{ar}{R}_{rnt}}-\frac{2B_{lmw}}{r\cdot R_{lnt}}\\ {D^{\prime }}_{13}=\frac{k_{al}{k}_{PT}}{R_{al}{R}_{lnt}}-\frac{k_{al}(k_{PR}-k_{DR})}{2d\cdot R_{al}{R}_{lnt}}-\frac{k_{ar}{k}_{br}}{r\cdot R_{ar}{R}_{rnt}}+\frac{B_{rmw}}{r\cdot R_{rnt}}\\ {D^{\prime }}_{14}=\frac{k_{al}{k}_{PT}}{R_{al}{R}_{lnt}}-\frac{k_{ar}(k_{PR}-k_{DR})}{2d\cdot R_{ar}{R}_{rnt}}+\frac{B_{lmw}}{r\cdot R_{lnt}}-\frac{k_{al}{k}_{bl}}{r\cdot R_{al}{R}_{lnt}}\\ {D^{\prime }}_{15}=m+\frac{2k_{ar}{k}_{DT}}{R_{ar}{R}_{rnt}}+\frac{2I_{lmw}}{r\cdot R_{lnt}},{D^{\prime }}_{16}=\frac{k_{al}{k}_{PT}}{R_{al}{R}_{lnt}}+\frac{I_{rmw}}{r\cdot R_{rnt}}\\ {D^{\prime }}_{17}=\frac{k_{al}{k}_{PT}}{R_{al}{R}_{lnt}}+\frac{I_{lmw}}{r\cdot R_{lnt}},{D^{\prime }}_{18}=0,{D^{\prime }}_{19}=0\end{array}$ (11)

From certain of the above equations, the following dynamic model of the mobile robot is obtained: $\left[\begin{array}{c}{\dot{x}}_P\\ {\dot{y}}_P\\ \dot{\theta }\\ \dot{u}\\ \dot{w}\end{array}\right]=\left[\begin{array}{cccc}u\cos \theta -aw\sin \theta & 0& 0& {\bar{\delta }}_x\\ u\sin \theta +aw\cos \theta & 0& 0& {\bar{\delta }}_y\\ w& 0& 0& 0\\ \frac{{D^{\prime }}_{12}}{{D^{\prime }}_{15}}u+\frac{mb}{{D^{\prime }}_{15}}{w}^2& \frac{{D^{\prime }}_{11}}{{D^{\prime }}_{15}}& 0& {\bar{\delta }}_u\\ \frac{{B^{\prime }}_{15}}{{B^{\prime }}_{11}}w+\frac{mb}{{B^{\prime }}_{11}}uw& 0& \frac{{B^{\prime }}_{14}}{{B^{\prime }}_{11}}& {\bar{\delta }}_w\end{array}\right]\left[\begin{array}{c}1\\ u_{ref}\\ w_{ref}\\ 1\end{array}\right]$ (13) $\begin{array}{c}{\bar{\delta }}_x=-0.5*({\dot{\eta }}_r+{\dot{\eta }}_l)\sin \theta \\ {\bar{\delta }}_y=0.5*({\dot{\eta }}_r+{\dot{\eta }}_l)\cos \theta \\ {\bar{\delta }}_w=\frac{{B^{\prime }}_{16}}{{B^{\prime }}_{11}}{\dot{\xi }}_r+\frac{{B^{\prime }}_{17}}{{B^{\prime }}_{11}}{\dot{\xi }}_l+(\frac{{B^{\prime }}_{19}}{{B^{\prime }}_{11}}-\frac{1}{2}mb){\ddot{\xi }}_r+\frac{{B^{\prime }}_{20}}{{B^{\prime }}_{11}}{\ddot{\xi }}_l\\ \begin{array}{cc}& \end{array}+\frac{e{F}_{e{y}^{\prime }}}{{B^{\prime }}_{11}}+\frac{c{F}_{c{y}^{\prime }}}{{B^{\prime }}_{11}}+\frac{{\tau }_e}{{B^{\prime }}_{11}}\\ {\bar{\delta }}_u=\frac{1}{{D^{\prime }}_{15}}({D^{\prime }}_{13}{\dot{\xi }}_r+{D^{\prime }}_{14}{\dot{\xi }}_l+{D^{\prime }}_{16}{\ddot{\xi }}_r+{D^{\prime }}_{17}{\ddot{\xi }}_l)\\ \begin{array}{cc}& \end{array}+\frac{1}{{D^{\prime }}_{15}}(F_{e{x}^{\prime }}+F_{c{x}^{\prime }}+\frac{1}{2}m{\dot{\eta }}_r+\frac{1}{2}m{\dot{\eta }}_l)\end{array}$

From equations (13), we can propose that the robotic acceleration variables u̇ and $\dot{\omega }$ do not depend on the position states x,y,and □ for robot; so the dynamics equation of this robot can be written as follows: $\left[\begin{array}{c}\dot{u}\\ \dot{w}\end{array}\right]=\left[\begin{array}{c}\frac{{D^{\prime }}_{12}}{{D^{\prime }}_{15}}u+\frac{mb}{{D^{\prime }}_{15}}{w}^2+\frac{{D^{\prime }}_{11}}{{D^{\prime }}_{15}}{u}_{ref}\\ \frac{{B^{\prime }}_{15}}{{B^{\prime }}_{11}}w+\frac{mb}{{B^{\prime }}_{11}}uw+\frac{{B^{\prime }}_{14}}{{B^{\prime }}_{11}}{w}_{ref}\end{array}\right]+\left[\begin{array}{c}{\bar{\delta }}_u\\ {\bar{\delta }}_w\end{array}\right]$ (14)

By re-arranging equations (14), we can obtain the following equations: $\left[\begin{array}{c}{u}_{ref}\\ w_{ref}\end{array}\right]=\left[\begin{array}{cccccc}\dot{u}& 0& u& 0& w^2& 0\\ 0& \dot{w}& 0& w& 0& uw\end{array}\right]\left[\begin{array}{c}{A}_1\\ A_2\\ A_3\\ A_4\\ A_5\\ A_6\end{array}\right]$ (15)

where: $\begin{array}{l}{A}_1=\frac{{D^{\prime }}_{15}}{{D^{\prime }}_{11}},A_2=\frac{{B^{\prime }}_{11}}{{B^{\prime }}_{14}},A_3=-\frac{{D^{\prime }}_{12}}{{D^{\prime }}_{11}},A_4=-\frac{{B^{\prime }}_{15}}{{B^{\prime }}_{14}}\\ A_5=-\frac{mb}{{D^{\prime }}_{11}},A_6=-\frac{mb}{{B^{\prime }}_{14}}\end{array}$

2.2 State description of MRS

In this section, the control laws developed so far are extended to the case of unicycle-type robot formations. The kinematic model of each agent or robot R_i(i = 1,…,n) is given by: $\left[\begin{array}{c}{\dot{x}}_{iP}\\ {\dot{y}}_{iP}\\ {\dot{\theta }}_i\\ {\dot{u}}_i\\ {\dot{w}}_i\end{array}\right]=\left[\begin{array}{cccc}{u}_i\cos {\theta }_i-a_i{w}_i\sin {\theta }_i& 0& 0& {\bar{\delta }}_{x_i}\\ u_i\sin {\theta }_i+a_i{w}_i\cos {\theta }_i& 0& 0& {\bar{\delta }}_{y_i}\\ w_i& 0& 0& 0\\ \frac{{D^{\prime }}_{12}}{{D^{\prime }}_{15}}{u}_i+\frac{m_i{b}_i}{{D^{\prime }}_{15}}{w}_i{}^2& \frac{{D^{\prime }}_{11}}{{D^{\prime }}_{15}}& 0& {\bar{\delta }}_{u_i}\\ \frac{{B^{\prime }}_{15}}{{B^{\prime }}_{11}}{w}_i+\frac{m_i{b}_i}{{B^{\prime }}_{11}}{u}_i{w}_i& 0& \frac{{B^{\prime }}_{14}}{{B^{\prime }}_{11}}& {\bar{\delta }}_{w_i}\end{array}\right]\left[\begin{array}{c}1\\ u_{ref}^{(i)}\\ w_{ref}^{(i)}\\ 1\end{array}\right]$ (16)

where u _ref ⁽ⁱ⁾ and w _ref ⁽ⁱ⁾ are, respectively, the instantaneous reference linear velocity and the reference angular velocity of the robot R_i, while u̇_i and ẇ_i are the instantaneous translation and rotational accelerations of the i-th robot. The state of the i-th mobile robot is then described by X_i = [x_ip y_ip θ _i u_i w_i] ^T , where i = 0 represents the virtual leader and i = 1,2,…, n represent other robots.

3. Artificial potential field functions

In this paper, a collision free trajectory of a multi-robot system under kinodynamic constraints in a fixed and bounded workspace is proposed. It is assumed that the robots perceive the whole environment. An adaptive formation control law for non-holonomic dynamic robots based on the artificial potential function method in the presence of lateral slip and parametric uncertainties is presented so that the robots can move safely towards their target points while maintaining the desired regular geometric formation.

3.1 Attraction to the target

In the paper, we can assume that a regular geometry formation centre position is assigned to the virtual leader while it arbitrarily moves on a collision-free trajectory. According to the formation geometry property, we can firstly define some instantaneous virtual targets related to the virtual leader which can be assigned to every robot. According to these virtual targets, they can be regarded as an attraction area of the artificial potential field which these targets locate in order to build an ordered formation in the desired position (the closest virtual target to themselves). The shape of the attraction area is a circle with a radius r_target. The value of the artificial potential in the attraction area is constant and grows outside the circle.

The attraction influence tends to pull the robot towards the target position area, while a repulsion force tends to push the robot away from the obstacles.

In a 2D planar, we consider an attractive potential function: $V_{att}^{(i)}(r_i)=\{\begin{array}{l}0\begin{array}{ccc}& & \begin{array}{cc}\begin{array}{cc}& \end{array}& \end{array}0<\Vert s_i\Vert \le r_{t\arg et}\end{array}\\ {\alpha }_{att}{(\Vert s_i\Vert -r_{t\arg et})}^{\frac{3}{2}}\begin{array}{cc}& \Vert s_i\Vert \ge r_{t\arg et}\end{array}\end{array}$ (17)

where α _att is the factor that determines the slope of the artificial potential field outside the attraction area, ‖s_i‖ is the distance between the i-th robot and the centre of the attraction area of its virtual target. The attraction force of i-th robot's virtual target is as follows: $f_{att}^{(i)}(r_i)=\{\begin{array}{l}0\begin{array}{ccc}& & \begin{array}{cc}\begin{array}{cc}& \end{array}& \end{array}0<\Vert s_i\Vert \le r_{t\arg et}\end{array}\\ \frac{3}{2}{\alpha }_{att}{(\Vert s_i\Vert -r_{t\arg et})}^{\frac{1}{2}}\begin{array}{cc}& \Vert s_i\Vert \ge r_{t\arg et}\end{array}\end{array}$ (18)

3.2 A repulsion potential field for obstacles

In order to make the mobile robots avoid all of the stationary obstacles intersecting their paths, we need to construct the obstacle avoidance functions which act close to the obstacles and decay as the distances which are measured by a laser range finder or vision system grow. It is assumed that there are q solid obstacles within the workspace and that the l-th obstacle is circular, with a centre position (o_l1,o_l2), and that r_ol is the radius of the forbidden area that surrounds l-th obstacle, q_ol is the radius of the area around the l-th obstacle where repulsion acts, $\Vert g_{il}\Vert$ is the distance between the i-th robot and the centre of the l-th forbidden area and α _obs determines the slope of the artificial potential field around the forbidden area.

When these robots are near to the obstacles, the artificial potential field that acts on the robots is proposed as follows [14]: $V_{obs}^{(i)}=\{\begin{array}{c}\begin{array}{cc}& \end{array}+\infty \begin{array}{cc}& \end{array}\begin{array}{ccc}& & \begin{array}{cc}& \end{array}0<\Vert g_{il}\Vert <r_{ol}\end{array}\\ \frac{{\alpha }_{obs}}{\Vert g_{il}\Vert -r_{ol}}+\Vert g_{il}\Vert \frac{{\alpha }_{obs}}{{(q_{ol}-r_{ol})}^2}\begin{array}{cc}& r_{ol}\le \Vert g_{il}\Vert <q_{ol}\end{array}\\ {\alpha }_{obs}\left(\frac{2q_{ol}-r_{ol}}{(q_{ol}-r_{ol})}\right)\begin{array}{ccc}& & \Vert g_{il}\Vert \ge q_{ol}\end{array}\end{array}$ (19)

where the inequality q_ol > r_ol must always be satisfied.

The artificial potential field between the robots and obstacles deduces the following repulsion force in order to avoid obstacles: $f_{obs}^{(i)}=\{\begin{array}{c}\begin{array}{cc}& \end{array}+\infty \begin{array}{cc}& \end{array}\begin{array}{ccc}& & \begin{array}{cc}& \end{array}0<\Vert g_{il}\Vert <r_{ol}\end{array}\\ -\frac{{\alpha }_{obs}}{{(\Vert g_{il}\Vert -r_{ol})}^2}+\frac{{\alpha }_{obs}}{{(q_{ol}-r_{ol})}^2}\begin{array}{cc}& r_{ol}\le \Vert g_{il}\Vert <q_{ol}\end{array}\\ 0\begin{array}{ccc}\begin{array}{cc}& \end{array}& & \end{array}\begin{array}{ccc}& & \Vert g_{il}\Vert \ge q_{ol}\end{array}\end{array}$ (20)

3.3 Interactions between robots

An artificial potential function is derived to control the interconnecting distance between a pair of mobile robots, and steer them towards their desired distances without collision. The formation distance error q_ij is defined as: $q_{ij}=\Vert {\overrightarrow{r}}_j-{\overrightarrow{r}}_i\Vert -d_{ij}$ (21)

where d_ij is the desired distance between the two robots and the i-th robot is denoted with the position vector ${\overrightarrow{r}}_i=(x_i,y_i)$ .

In the formation control of MRS, the artificial potential fields around the robots maintain the MRS formation with a desired regular geometry in order to avoid collisions among them. The strength of these forces -which are caused by the artificial potential fields among the robots – depends upon the distance between the formation members. If these robots are too close, these robots will be mutually repulsed and if these robots are too far, they will be mutually attracted. When the distance between the robots exceeds a specified limit value, the attraction forces will become zero. In order to present the attraction potential and the repulsion potential, the potential energy between the i-th and the j-th robots is constructed as follows: $V_{\mathrm{int}er}^{(i,j)}(r_i,r_j)={\alpha }_{\mathrm{int}er}^{(i,j)}\ln (\cosh (q_{ij}))$ (22)

where α_inter ^i,j determines a slope of the artificial potential around the robot: $f_{\mathrm{int}er}^{(i,j)}(r_i,r_j)=\nabla V_{\mathrm{int}er}^{(i,j)}(r_i,r_j)={\alpha }_{\mathrm{int}er}^{(i,j)}\tanh (q_{ij})\frac{r_j-r_i}{\Vert r_j-r_i\Vert }$ (23)

4. Formation control law based on artificial potential field functions

4.1 Formation control law

The control law based on artificial potential fields for the i-th robot combines the repulsion to those obstacles, the attraction to the virtual target's attraction areas belonging to i-th robot and the interaction between the robots, as proposed by the following equation: $U_{total}^{(i)}(r_i)=V_{att}^{(i)}(r_i)+{\displaystyle \sum _{j=1,j\ne i}^n{V}_{\mathrm{int}er}^{(i,j)}}(r_i,r_j)+{\displaystyle \sum _{l=1}^q{V}_{obs}^{(i,l)}(r_i)}$ (24)

The control of non-holonomic mobile robot information is based on the incorporation of all of the potential forces. The synthesized force for the i-th robot, F _total ⁽ⁱ⁾ , is given: $F_{total}^{(i)}(r_i)={\nabla }_{r_i}{U}_{total}^{(i)}$ (25)

where ∇ _{r i} U_total⁽ⁱ⁾ includes the potential energy from the destination target, the obstacles spread in the environment and the interactions between robots.

The main objective of the formation control law is to assign a group of non-holonomic mobile robots to form a regular geometric shape. The descending gradient vector field is used as a reference for the heading direction of the robot; the projection of $-{\nabla }_{x_i}{U}_{total}^{(i)}(r_i,r_j)$ and $-{\nabla }_{y_i}{U}_{total}^{(i)}(r_i,r_j)$ onto the current heading direction of the robot is used to control its translation acceleration for the robot.

In order to arrive at a desired target point and form a desired formation shape in an unknown environment, the formation control law tries to reduce the orientation error and asymptotically converge on zero or a minimum value while maintaining the desired formation shape and arriving at the target position. Because the angle of the synthesized potential force can act as a trajectory of the mobile robot, the orientation error of the i-th robot evolves over time, as follows: $\{\begin{array}{c}{e}_i={\theta }_i-{\phi }_i\\ {\dot{e}}_i={\dot{\theta }}_i-{\dot{\phi }}_i\\ {\ddot{e}}_i={\ddot{\theta }}_i-{\ddot{\phi }}_i\end{array}$ (26)

Theorem 1: Consider a multi-robot system whose motion is governed by the Ordinary Differential Equation (ODE) described in Equation (13). If we assume that the synthesized force defined by the target's attraction force, the obstacle's repulsion force and the interactions between the robots is applied to each robot and that the angle of the synthesized force in the global coordinates [16], ${\phi }_i=\mathrm{atan}2(-\frac{\partial U_{total}^{(i)}}{\partial y_i},-\frac{\partial U_{total}^{(i)}}{\partial x_i}),$

is regarded as the desired heading of each robot's movement.

In order to stabilize a prescribed formation, facilitate the manoeuvres of the robots within a constrained environment and asymptotically reach the target configuration while maintaining a desired formation, the decentralized control law is established as follows: $\begin{array}{l}{u}_{ref}^{(i)}=\frac{1}{u_i{k}_{u_i}}\frac{{D^{\prime }}_{15}}{{D^{\prime }}_{11}}\left[-f_x^{(i)}{u}_i\cos {\theta }_i+f_x^{(i)}{a}_i{w}_i\sin {\theta }_i\right]\\ \begin{array}{cc}& \end{array}+\frac{1}{u_i{k}_{u_i}}\frac{{D^{\prime }}_{15}}{{D^{\prime }}_{11}}\left[-f_y^{(i)}{u}_i\sin {\theta }_i-f_y^{(i)}{a}_i{w}_i\cos {\theta }_i\right]\\ \begin{array}{cc}& -\end{array}\left[2\frac{{D^{\prime }}_{12}}{{D^{\prime }}_{11}}{u}_i+\frac{m_i{b}_i}{{D^{\prime }}_{11}}{w}_i^2+\frac{{D^{\prime }}_{15}}{{D^{\prime }}_{11}}{\bar{\delta }}_{u_i}\right]\end{array}$ (27) $\begin{array}{l}{w}_{ref}^{(i)}=-\frac{{B^{\prime }}_{11}}{{B^{\prime }}_{14}}{k}_{p_i}{e}_i-2\frac{{B^{\prime }}_{15}}{{B^{\prime }}_{14}}{k}_{d_i}{\dot{e}}_i-\frac{{B^{\prime }}_{15}}{{B^{\prime }}_{14}}{k}_{d_i}{\dot{\phi }}_i-\frac{m_i{b}_i}{{B^{\prime }}_{14}}{u}_i{\dot{e}}_i\\ \begin{array}{cc}& \end{array}-\frac{m_i{b}_i}{{B^{\prime }}_{14}}{u}_i{\dot{\phi }}_i+\frac{{B^{\prime }}_{11}}{{B^{\prime }}_{14}}{\ddot{\phi }}_i-\frac{{B^{\prime }}_{11}}{{B^{\prime }}_{14}}{\bar{\delta }}_{w_i}\end{array}$ (28)

where k_di, k_pi, k_ui are positive gains and f_x⁽ⁱ⁾, f_y⁽ⁱ⁾ are the components of the synthesized artificial potential force on the i-th robot in the x and y directions, respectively.

Proof: Consider the following Lyapunov candidate [15]: $W(X)=\frac{1}{2}{\displaystyle \sum _{i=1}^n{U}_{total}^{(i)}}(r_i)+\frac{1}{2}{\displaystyle \sum _{i=1}^n{k}_{p_i}}{e}_i^2+\frac{1}{2}{\displaystyle \sum _{i=1}^n{\dot{e}}_i^2}+\frac{1}{2}{\displaystyle \sum _{i=1}^n{k}_{u_i}{u}_i^2}$ (29)

Cleary, the candidate function satisfies the prerequisites of a Lyapunov function in that: (1) if, and only if, we defined X_i: (s_i,e_i,ė_i,u_i,w_i) = 0 as an equilibrium point of the MRS formation control problem, we can obtain W(X) = 0, and (2) W(X) is a locally positive and continuously differentiable function on the domain $D(W)=\left\{X\in {ℝ}^{5\times n}:V_{obs}^{(i)}>0\right\}$ . The time derivative of our Lyapunov function W along a particular trajectory of the system Equation (13) is computed as: $\dot{W}(X)={\displaystyle \sum _{i=1}^n{({\nabla }_{r_i}{U}_{total}^{(i)})}^T}{\dot{r}}_i+k_{p_i}{e}_i{\dot{e}}_i+{\dot{e}}_i{\ddot{e}}_i+k_{u_i}{u}_i{\dot{u}}_i$ (30) $\begin{array}{l}{\left({\nabla }_{r_i}{U}_{total}^{(i)}\right)}^T{\dot{r}}_i=\left(\begin{array}{cc}{({\nabla }_{r_i}{U}_{total}^{(i)})}_x& {({\nabla }_{r_i}{U}_{total}^{(i)})}_y\end{array}\right)\left(\begin{array}{c}{u}_i\cos {\theta }_i-a_i{w}_i\sin {\theta }_i\\ u_i\sin {\theta }_i+a_i{w}_i\cos {\theta }_i\end{array}\right)\\ \begin{array}{cccc}& & & \end{array}=(\begin{array}{cc}{f}_x^{(i)}& f_y^{(i)}\end{array})\left(\begin{array}{c}{u}_i\cos {\theta }_i-a_i{w}_i\sin {\theta }_i\\ u_i\sin {\theta }_i+a_i{w}_i\cos {\theta }_i\end{array}\right)\end{array}$

Inserting the formation control law equations (27), (28) and the system kinematics (13) into Equation (30), we get: $\begin{array}{l}\dot{W}(X)={\displaystyle \sum _{i=1}^n[(\begin{array}{cc}{f}_x^{(i)}& f_y^{(i)}\end{array})\left(\begin{array}{c}{u}_i\cos {\theta }_i-a_i{w}_i\sin {\theta }_i\\ u_i\sin {\theta }_i+a_i{w}_i\cos {\theta }_i\end{array}\right)}\\ \begin{array}{ccc}& & +k_{p_i}{e}_i\end{array}{\dot{e}}_i+{\dot{e}}_i({\dot{w}}_i-{\ddot{\phi }}_i)+k_{u_i}{u}_i{\dot{u}}_i\begin{array}{c}\\ \end{array}]\\ ={\displaystyle \sum _{i=1}^n[\begin{array}{c}\\ \end{array}(-f_x^{(i)}{a}_i{w}_i\sin {\theta }_i+f_x^{(i)}{u}_i\cos {\theta }_i+f_y^{(i)}{u}_i\sin {\theta }_i+f_y^{(i)}{a}_i{w}_i\cos {\theta }_i)}\\ +k_{p_i}{e}_i{\dot{e}}_i+{\dot{e}}_i\left(\frac{{B^{\prime }}_{15}}{{B^{\prime }}_{11}}{w}_i+\frac{m_i{b}_i}{{B^{\prime }}_{11}}{u}_i{w}_i+\frac{{B^{\prime }}_{14}}{{B^{\prime }}_{11}}{w}_{ref}^{(i)}+{\bar{\delta }}_{w_i}-{\ddot{\phi }}_i\right)\\ +k_{u_i}{u}_i\left(\frac{{D^{\prime }}_{12}}{{D^{\prime }}_{15}}{u}_i+\frac{m_i{b}_i}{{D^{\prime }}_{15}}{w}_i^2+\frac{{D^{\prime }}_{11}}{{D^{\prime }}_{15}}{u}_{ref}^{(i)}+{\bar{\delta }}_{u_i}\right)\begin{array}{c}\\ \end{array}]\\ ={\displaystyle \sum _{i=1}^n(-k_{d_i}{\dot{e}}_i^2-k_{u_i}{u}_i^2)\le 0}\end{array}$ (31)

According to the definition of the Lyapunov candidate function W(X) ≥ 0 and the result whereby W(X) ≤ 0, it is implied that the function W(X) is a legitimate Lyapunov function guaranteeing the stability and asymptotic approach of the desired formation [14].

5. Multi-target optimal assignment for the formation problem

The optimal assignment (distribution) problem is encountered in the formation of MRSs. The following two assumptions should be satisfied when n robots are appointed to reach m target positions which are predefined by the desired formation geometry shape:

If and only if each robot reach a target point;

If and only if each target point be occupied by one robot.

where C_ij is the resource consumed when the robot i reaches the j-th target point.

There are various ways to solve the above assignment problems. Certainly it can be formulated as a linear program and the simplex method can be used to solve it. However, sometimes the simplex method is inefficient for the formation of MRS assignment problems. The improved Hungarian algorithm has been used to solve these formation assignment problems and is summarized as follows: Step 1. Determine the cost table for the given problem;

If the total of the target points m is equal to the total of robots n, go to Step 3;

Otherwise, go to Step 2.

Step 2. If n < m, add m-n dummy robots so that the cost table becomes a square matrix. The cost entries of the dummy robots are always zero. If n > m, use the extended contract net protocol [19] to choose m robots from n candidates to reach target point m in formation.

Step 3. Locate the smallest element in each row of the cost matrix and then subtract the same from each element of the row.

Step 4. In the reduced matrix obtained by Step 3, locate the smallest element of each column and then subtract the same from each element of that column. Each column and row now have at least one zero.

Step 5. In the new matrix obtained from Step 4, search the optimal assignment as follows:

• a)

  Test every row of the matrix until a row with a minimum of zero is found. Mark this row with a particular symbol (for example, a small rectangle symbol‘□’) and cross off all the other zeros in its column with symbols (X). Continue with this process until all of the rows have been inspected.

• b)

  Repeat the procedure for each column of the reduced matrix.

• c)

  If a row and/or a column have two or more zeros and one cannot be chosen by inspection, then assign arbitrarily any one of these zeros and cross off all the other zeros of that row/column.

• d)

  Repeat (a), (b) and (c) until the chain of assignment (□) or cross (X) ends.

Step 6. If the number of the assignment (□) is equal to n, an optimum solution is reached. If the number of assignments is less than n, go to the next step 7.

Step 7. Draw the minimum number of horizontal and/or vertical lines to cover all the zeros of the above reduced matrix.

Step 8. Obtain the new modified cost matrix as follows:

• a)

Research the smallest element of the modified matrix not covered by any of the lines.

b)

Subtract this element from all of the uncovered elements and add the same to all of the elements lying at the intersection of any two lines.

Step 9. Go to Step 6 and repeat the procedure until an optimum solution is reached.

6. Evaluation method of formation shape control

In order to quantitatively evaluate the formation shape control and capture the shape distortion in our paper, the concept of Procrustes shape distance [22] is used. The Procrustes shape distance is defined based on the following definitions.

Definition 1. The Configuration is the set of landmarks on a particular object. The vector representation X is the n × k matrix of the Cartesian coordinates of n landmarks in k dimensions which represent the configuration matrix.

$X=\left[\begin{array}{cccc}{x}_{11}& x_{12}& \cdots & x_{1k}\\ x_{21}& x_{22}& \cdots & x_{2k}\\ ⋮& ⋮& ⋮& ⋮\\ x_{n1}& x_{n2}& \cdots & x_{nk}\end{array}\right]$

Definition 2. Shape is all of the geometrical information that remains when the location, scale and rotational effects are filtered out from an object. Let x₁,…,x_n and y₁,…,y_n be two (not totally coincident) labelled sets of n points in ℛ^k. We regard these two configurations as having the same shape if y_i = cRx_i + b,i = 1,2,…,n for some scale factor c > 0, a rotation matrix R(k × k),and a vector b(n × 1). In order to measure distances between shapes (or distortion shapes) we require a metric – for example, Procrustes distance as the metric.

For a description of shape in coordinate terms, consider a configuration of n points in ℛ^k, represented by a n × k matrix X. Location and scale effects are easy to eliminate directly. As such, we may recall that the standard Helmert matrix has its first row of elements equal to $n^{-\frac{1}{2}}$ and the remaining rows are orthogonal to this first row. Let this Helmert matrix minus its first row be a Helmert sub-matrix which is satisfied as follows: ${(H)}_j=(d_j,\mathrm{\dots },d_j,-j{d}_j,0,\mathrm{\dots },0),d_j={\left\{j(j+1)\right\}}^{-\frac{1}{2}},j=2,\mathrm{\dots },n$

where dj is repeated j times. The rows of the Helmert matrix H are orthonormal and H1 _n = 0, 1 _n represents a n-vector of one. As such, Y = HX(n -1 × k) is invariant under the location shifts of the configuration.

The Helmert sub-matrix H can be used to remove the location of a configuration X. The full Helmert matrix H^F is a (n × n) orthogonal matrix with all of the elements in its first row equal to $\frac{1}{\sqrt{n}}$ and all of the remaining rows orthogonal to the first row, such that: $H^F=\left[\begin{array}{ccc}\frac{1}{\sqrt{3}}& \frac{1}{\sqrt{3}}& \frac{1}{\sqrt{3}}\\ -\frac{1}{\sqrt{2}}& \frac{1}{\sqrt{2}}& 0\\ -\frac{1}{\sqrt{6}}& -\frac{1}{\sqrt{6}}& \frac{2}{\sqrt{6}}\end{array}\right],H=\left[\begin{array}{ccc}-\frac{1}{\sqrt{2}}& \frac{1}{\sqrt{2}}& 0\\ -\frac{1}{\sqrt{6}}& -\frac{1}{\sqrt{6}}& \frac{2}{\sqrt{6}}\end{array}\right]$

Scale effects can be eliminated by centring the configuration and then dividing by size, as follows: $Z=\frac{HX}{\Vert HX\Vert }$

Accordingly, Z can be regarded as a point on the unit sphere in ℛ^(n-1)k, which can be described as the pre-shape matrix.

One route to obtaining shape space is to first remove the effect of location and scale, giving “pre-shape”.

If we remove translation from the original configuration, then the resulting landmarks are called Helmertized. Filtering scale from these Helmertized landmarks yields pre-shape whereas eliminating rotation from them should create size-and-shape. Again, removing rotation from pre-shape or removing scale from size-and-shape should result in shape. The size-and-shape configuration [X] _s can be represented as: ${[X]}_S=\left\{X_H\mathrm{\Gamma }:\mathrm{\Gamma }\in SO(k)\right\}$

where X_H = HX are the Helmertized coordinates and γ is an (k × k) rotation matrix in the special orthogonal group SO(k). Finally, if size is removed from the size-and-shape of a configuration X by re-scaling to unit centroid size, then the resulting configuration describes the shape of X. The shape can be described by the set [X] given by: $[X]=\left\{Z\mathrm{\Gamma }:\mathrm{\Gamma }\in SO(k)\right\}$

The Procrustes distance dp is obtained by matching the Helmert coordinates X_Hr and X_Hd of X_r and X_d as closely as possible over rotations. $d_p(X_r,X_d)=\underset{\mathrm{\Gamma }\in SO(k)}{\inf }\Vert X_{Hr}-X_{Hd}\mathrm{\Gamma }\Vert$

Meanwhile, this Procrustes distance can be used to represent shape distortion in our formation control task, where X_r is the configuration matrix for three WMRs in the actual formation (in order to compare with reference [20], we proposed the simplest rectangular formation with three robots) and X_d is the configuration matrix for three WMRs in the desired formation. Both X_r and X_d are 3×3 matrices.

7. Experiments and results

Experiments to test the control of robot formations were both performed using three erratic mobile robots and a Pioneer mobile robot and a simulation platform, where the Pioneer mobile robot was equipped with a laser range finder and CCD with Pan-Tilt-Zoom (PTZ) control. Every Erratic robot equipped with both an omni-directional camera system mounted upon it and a URG-04LX laser range finder with an on-board computer for image processing (see Fig. 2).

Figure 2.

Four Mobile robots and an omni-directional camera and laser finder.

Simulations were performed with up to ten robots. Four regular geometry formation shapes (for example, triangle, circle, line and diamond formations) were demonstrated in the experiment.

7.1 Experiment on a real multi-robot platform

An off-board computer is used to monitor and start the formation of the multi-robot system through a wireless Ethernet and radio modem.

During formation forming and maintenance, the identification of the role was performed using a colour marker.

In Figure 3, a sequence of video snapshots of the experiments in HIT Multi-Agent robotic research centre are presented, demonstrating formation generation and formation shrinking with three Erratics robots with an omni-directional camera system and a Pioneer mobile robot (playing the role of the static virtual leader).

Figure 3.

Experiments on formation generation in MRS: (a) Initial set-up, (b) Formation generation in process, (c) Formation shrinkage, (d) Formation generation

In Fig.4, the virtual leader moves with a random strategy, the three Erratics robots' initial formation is a line in the formation generation process with the three robots following the leader robot (Pioneer robot) and generating a triangular formation closing on the leader.

Figure 4.

Formation generation in MRS (the leader is moving): (a) Initial line formation set-up, (b) Fan formation generation in process, (c) Formation shrinkage, (d) Formation generation

Formation generation, formation maintenance and formation change at the MRS simulation platform constitute the three scenarios considered in this paper. In the first scenario (formation generation) a number of robots will form a regular polygon – for example, for triangle, circle, line and diamond formations, etc. The whole robot team uses artificial potential fields to deal with obstacle avoidance and goes to the target point in the complex environment. This is related with the second scenario – formation maintenance – in which the virtual leader can move arbitrarily in this environment while avoiding obstacles and other robots. In the third scenario, the MRS can finish the formation task by changing among the above four formations.

7.2 Formation generation on the simulation platform

In this section, different robot teams will be considered. The teams vary in size from three up to twelve, in which situation a robot can be homogeneous or heterogeneous with different capabilities. We consider three to be the minimum number of robots to form a regular polygon formation around the virtual leader. Meanwhile, we simulate a maximum number of 12 robots, since a higher number of robots will demonstrate similar behaviours. Equations (28–29) are used as the control model governing each robot. In Fig. 5 can be seen the trajectories of different robot teams forming a regular polygon around the virtual leader. The initial position of each robot is random. The black areas represent different sized obstacles. In Fig. 5, different robot teams end up in a regular polygon formation around the virtual leader. Just as with Fig. 5, if there is enough space for the robots to form a regular polygon – even though there are some obstacle on the path from the initial position to the desired target position for a robot – the robots always converge into a stable situation using the artificial potential field method.

Figure 5.

Trajectories of different simulated robot teams forming a regular polygon around a virtual leader with (a) A triangular formation with 9 robots, (b) A circular formation with 9 robots, (c) A triangular formation with 3 robots, (d) A circular formation with three robots, (e) A diamond formation with 8 robots, (f) A diamond formation with twelve robots

7.3 Formation maintenance on the simulation platform

In this section, we simulate different robot teams in an environment with obstacles aiming to maintain the desired regular polygon when the virtual leader moves arbitrarily in the environment while avoiding obstacles and other robots.

Fig. 6 and Fig. 7 present the situation where the virtual leader moves following the small rectangle: how do the robot teams maintain their regular formation? At the initial instance, all of the robots will be affected by the attractive force from the desired target point and the repulsive force from any obstacle nearby while another robot will plan its own path free from obstacles. When they reach their own desired target point, they have to dynamically maintain their regular geometry formation while following the virtual leader's movement. From Fig. 6 and Fig. 7, it can be seen that when the virtual leader's trajectory is a small rectangle, all of the robots travel in the same, almost rectangular trajectory.

Figure 6.

Simulation experiment on Line formation keeping in a seven robot team's formation with obstacles in the environment at different moment

Figure 7.

Simulation experiment on a circular formation maintaining an 8 robot team formation with obstacles in the environment at different moments

7.4 Formation dynamic changing among four regular geometry formations

Four regular geometry formations can be mutually transformed, the transformation processes of which are as follows (Fig. 8): when the robot team forms a regular formation, if we want to change formation to another regular formation, we only have to press one button on the keyboard to finish the change runtime.

For example: “p” or “P” represents the transformation from the current formation to a diamond formation, “c” or “C” represents the transformation from the current formation to a circle formation, “l” or “L” represents the transformation from the current formation to a line formation, “t” or “T” represents the transformation from the current formation to a rectangle formation. The proposed algorithm will recalculate the new target position according to the new formation and the new position of the virtual leader and reassign the above targets to robots using the improved Hungarian algorithm. Once one robot obtains its target position, it will plan its path from its position to the target position without obstacles using APF and go to the target following the control law equation (Eq. 28–29).

Figure 8.

Formation dynamic changing among four regular geometry formations with 8 robots: (a) Initial line formation, (b) Change from line to circle formation, (c) Change from circle to rectangle formation, (d) Change from rectangle to diamond formation, (e) Change from diamond to line formation, (f) Change from line to rectangle formation

7.5 Formation shape distortion without slip

In Fig. 9, we observe the path of the three WMRs following a rectangular path (or equal to four L-shapes) while preserving a triangular formation. The blue triangular symbol denotes the following robot 2's position, the red circle symbol denotes the following robot 1's position, the green rectangle symbol represents the leader's position.

Figure 9.

Rectangle shape formation for three WMRs without wheel-slip

Fig. 10 shows the shape distortion for three WMR's rectangular path while trying to preserve the desired formation. As can be observed from [20], initially there was a shape distortion which the MRS system was trying to reduce before the leader took a sharp turn (about 90 degrees). When the leader first turned sharply, the shape distortion curve of the following robot 1 formed an initial peak (about 80∼115 s) due to robot 1's need to sharply change its own position from a position below the leader to a position left of the leader; when the leader completes the second sharp turn (about 180 degrees), the curve meets the second peak and the following robot 2 needs to sharply change its position from a position right of the leader to a position under the leader; when the virtual leader moves to a closed rectangular path, the curve has four peaks (every following robot is subject to two distortion peak values).

Figure 10.

Shape distortion for three WMRs without wheel-slip

In Fig. 11, we observe the path of the three WMRs following a rectangular path (or equal to four L-shapes) while preserving a triangular formation using the Hungarian algorithm to optimally assign a formation target position. From Fig. 12, we can observe that the following robots can change their role dynamically by the optimal assignment – it will move to the closest target while maintaining the desired formation shape so as to reduce the distortion shape. When the virtual leader turns sharply (at about 80 s, 170 s, 260 s and 350 s), the shape distortion for the three following robots is very small (less than 5 m). From this figure, it can be seen that the optimal assignment algorithm for the formation is very efficient for reducing the shape distortion.

Figure 11.

Three WMRs' formation with one virtual leader

Figure 12.

Shape distortion for three WMRs using optimal assignment without wheel-slip

7.6 Formation shape distortion for MRS subject to wheel-slip

In the simulation experiment, we show that when the MRS meets a large amount of slip, the stable formation becomes unstable. In this case, we can define the friction coefficient as 0.3 for this simulation. The triangular formation evolution and shape distortion during the rectangular path following can be seen in Fig. 13 and Fig. 14. Because of wheel-slip, the leader cannot take a sharp turn or form a rectangular shape and, as a result, the formation control becomes unstable quickly. From Fig. 14, we can observe that the distortion shape is very large (an average of more than 50) and meets four distortion peak values such that the formation is extremely unsteady and difficultly to maintain.

Figure 13.

Rectangle shape Formation for three WMRs subject to wheel-slip

Figure 14.

Shape distortion for three WMRs subject to wheel-slip

In order to solve wheel-slip problem, we used a dynamic WMR model subject to wheel-slip so as to design a slip-sensitive controller (detailed in Section 2 and Section 3). From Fig. 15 and Fig. 16, it is possible to compensate the slip effectively and stabilize the formation control, for the convergence of the formation to the desired shape and the shape distortion to zero. Although the path of the virtual leader is not a standard rectangle due to the wheel-slip, the three WMRs' distortion shape is very small and quickly converges to a low value (less than 5).

Figure 15.

Three WMRs' formation using the optimal assignment subject to slip

Figure 16.

Shape distortions of WMRs using the optimal assignment subject to slip

8. Conclusion

In this article, a strategy for the formation control of a non-holonomic multi-robot system based on artificial potential functions while subject to wheel-slip was demonstrated. A complete dynamic model of a unicycle-like mobile robot with wheel-slip and its linear parameterization was presented. In addition, the controller for generating, transforming and maintaining formation around a virtual leader by forcing the motion of the robots along the gradient field of the virtual potential function generated by the inter-individual distance requirements. The simulation results show that the controller constitutes a suitable method for MRS when subject to wheel-slip formation and navigation and obstacle avoidance, and that it is robust for tracking desirable formations with a regular geometry evolving under a dynamic environment. The controller for the formation shape using optimal assignment can solve the distortion shape in the presence of some wheel-slip and avoid the robots' frequently changing their positions and roles.

In the future, more complex scenarios such as environments with more obstacles and narrow corridors will be considered. An algorithm to make robots use different control behaviours under more realistic kinematic and dynamic constraints will be included.