General Path Planning Methodology for Leader-Follower Robot Formations

This paper describes a robust algorithm for mobile robot formations based on the Voronoi Fast Marching path planning method. This is based on the propagation of a wave throughout the model of the environment, the wave expanding faster as the wave's distance from obstacles increases. This method provides smooth and safe trajectories and its computational efficiency allows us to maintain a good response time. The proposed method is based on a local-minima-free planner; it is complete and has an O(n) complexity order where n is the number of cells of the map. Simulation results show that the proposed algorithm generates good trajectories.


Introduction
Formation control is currently a highly popular topic of research. Different approaches can be classified according to different criteria. Beard et al. [1] consider the different design approaches and classify them into three groups. First, leader-following, where one vehicle is designated as the leader and the others as followers. The leader posture (position and orientation) is determined by a trajectory to be tracked or by external control objectives (e.g., joysticked by a human) and the followers must track the leader following some prescribed geometry, possibly dynamically changing over time [2]. Second, behavioural, where the motion of each vehicle results from a weighted average of several behaviours, ultimately contributing to a desired group behaviour [3]. Third, virtual structure, where the entire formation is treated as a single structure whose desired motion is translated into the desired motion of each vehicle [4]. Another possible criterion takes into account the rigidity of the formation geometry: some authors specify the full geometry, e.g., the distances and bearings between the vehicles of the formation and control each vehicle to ensure that these are accurately achieved [5], requiring a coordination architecture to switch between geometries when required by environment characteristics (e.g., narrow passages, open spaces) [6]; others see the formation as a dynamic geometry structure that naturally becomes distorted in the presence of obstacles and/or environment geometry changes [7].
We are particularly interested in leader-following deformable formations, where the leader can be virtual and tracks a given trajectory, pulling the followers behind it according to nominal geometry specifications (e.g., desired inter-vehicle distances) that can change within a given range so as to accommodate environment conditions. Leonard and Fiorelli [8] introduced the concept of artificial potential fields between formation vehicles, some of them virtual leaders. The nominal intervehicle distance corresponds to minima of the potential functions, representing equilibrium points that balance inter-vehicle repulsion, vehicle-obstacle repulsion, and follower-leader attraction. MacArthur and Crane proposed a similar approach, but used virtual springdamper systems to "connect" the formation vehicles [9]. The drawback of such approaches is the well-known local minima problem of potential fields, which may lead to the formation being broken in the presence of non-convex shaped obstacles.
In a previous paper [10], we have introduced one possible solution to this problem, where the followers keep track of the most recent positions of the leader, instead of only its current position, in order to be dragged away from obstacle traps. However, this method has not been proven to work for all possible situations. We have also introduced a method to avoid potential fields local minima for a single vehicle, using the Voronoi Fast Marching method (VFM) and the Fast Marching Squared method ( ) [11], [12], [13].
In this paper we use the Fast Marching (FM) algorithm to control a leader-follower deformable formation, where the trajectory of the leader in an environment cluttered by obstacles is computed using the algorithm. At each iteration step, each follower attempts to reach a nominal sub-goal position related to the desired leader trajectory, but takes into account the positions of the other followers and the environment objects, both seen as obstacles. This influences the metrics used by these algorithms, effectively deforming the followers' trajectories. This way we ensure that non-convex obstacles do not break the formation and that the inter-vehicle distances are smoothly deformed while the formation moves from open areas to regions with obstacles, narrow corridors, and narrow passages.
The paper is organized as follows: in Section 2 the Fast Marching method and some of its variants are summarized. Its application to robot formations is described in Section 3. The results of applying the method to two different formation geometries in a diversity of simulation scenarios are presented in Section 4. The paper ends with conclusions and prospects for future work (Section 5).

Introduction to Fast Marching and Level Sets
The FM algorithm was introduced by J. Sethian in 1996 and is a numerical algorithm that approximates the viscosity solution of the Eikonal equation (1) which describes for instance light propagation in a nonhomogeneous medium. The level set of the solution represents the wave front advancing with a medium velocity , which is the inverse of the medium refraction index . Therefore, the Eikonal equation can be written as . The resulting function is a distance function, and if the medium velocity is constant, it can be seen as the Euclidean distance function to a set of starting points, usually the goal points. If the medium is non-homogeneous and the velocity is not constant, the function represents the distance function measured with the metrics or the arrival time of the wave front to point .
The FM method is used to solve the Eikonal equation and is very similar to the Dijkstra algorithm that finds the shortest paths on graphs, though it is applied to continuous media. Using a gradient descent of the distance function , one is able to extract a good approximation of the shortest path (geodesic) in various settings (Euclidean distance with constant and a weighted Riemannian manifold with varying ).

Intuitive Introduction to the Fast Marching Motion Planner
There are attractive potentials for obtaining a motion planner for mobile robots with desirable properties such as smoothness and safety. Some natural phenomena show a similar behaviour, e.g., electromagnetic waves. If there is an antenna in the goal point that emits an electromagnetic wave, then the robot can drive to the destination by tracing the waves back to the source. In general, the concept of electromagnetic waves is especially interesting, since the potential and its associated vector field have the desired trajectory properties, including smoothness and the absence of local minima.
Nevertheless, there are some problems. The most important problem that typically arises in mobile robotics is that optimal motion plans may bring robots too close to obstacles or people, which is not safe. This problem has been addressed by Latombe [14], and the resulting navigation function is called . The Voronoi method also tries to follow a maximum clearance map [15]. To obtain a safe path, it is necessary to add a component that repels the robot from obstacles, as proposed by Khatib [16]. In addition, this repulsive potential and its associated vector field should have good properties such as those of electrical fields. If we consider that the robot has an electrical charge of the same sign as obstacles, then the robot will be pushed away from obstacles. The properties of this electric field are very good: it is smooth and there are no singular points in the interest space ( ).
The main problem of attractive and repulsive potentials is how to mix the two fields together. This union between an attractive and a repulsive field has been the biggest problem of potential fields in path planning since Khatib [16]. The problem is that the sum, difference, or similar operations provoke the appearance of local minima.
In the proposed algorithm, this problem has been solved in the same way as it is in nature: the electromagnetic waves, as the light, have a propagation velocity that depends on the medium due to the refraction index. For example, a flint glass has a refraction index of 1.6, while in the air it is approximately 1. This refraction index of the medium is the quotient between the light velocity in the vacuum and its velocity in the medium, i.e., the slowness index of the front wave propagation in a medium. This idea of mixing the two potentials as in nature is the main contribution of the VFM and FM � methods [11][12][13].
According to Fermat's principle, the path taken by light from one point to another is the path that minimizes the time. In the case of a homogenous medium in which the light velocity is constant, the light follows a straight line. The set of points achieved in a fixed time is a circle. Consider two media with different refractive indexes, e.g., air and water. If a ray of light passes from the first medium to the second, the phenomenon known as refraction occurs. The beam appears to bend when the medium changes. The light seems to prefer to stay more in the medium with greater light velocity, as shown in Figure 2a. If there is a continuous change in the refraction index of the medium, the path obtained is a curve that moves away from areas with lower light velocity (higher refractive index), as shown in Figure 2b and Figure 1. This is the case of the well-known hot road surface mirage, in which "fake water" on the road is produced by the light rays bending due to a change of the refraction index in air with higher temperature (and lower refraction index) near the road surface.  Figures 2c and 2d show the behaviour of the attractive potential used alone. When the attractive potential uses the repulsive potential as refraction index, a unique field is obtained and its associated vector field is attracted to the goal point and repulsed from the obstacles, as shown in Figures 2e and 2f. This method inherits the properties of the electromagnetic field, i.e., � � is found if the refraction index is � � . Intuitively, the FM method gives the propagation of a front wave in nonhomogenous media. This repulsive potential can be obtained using the Extended Voronoi Transform (EVT) of the binary image of the map. The EVT computes the Euclidean Distance of the binary image. In a binary image, a pixel is referred to as background if it corresponds to a complete absence of obstacles and its value is one. The value of a pixel is zero if it corresponds to an obstacle or a wall. For a given distance metric, the EVT of an image produces a distance map of the same size. For each pixel in the image, the EVT assigns a number which is the distance to the nearest zero pixel of the image. For each pixel inside the objects in the binary image, the corresponding pixel in the distance map has a value equal to the minimum distance to the background. This is the solution adopted in the VFM method.
Clearly, the EVT is closely related to the Voronoi diagram. The Voronoi diagram concept is involved in many EVT approaches, either explicitly or implicitly.
Another possibility is to build the repulsive potential using FM. This is the solution adopted in the FM � method. In this case, a wave is propagated starting from the points representing the obstacles and walls. This wave propagation is achieved through the FM method. The result is a potential map, which can be interpreted as a velocity map (or slowness map) because it gives a clear idea of the permissible robot velocity at each environment cell. This potential map is represented in grey scale, where the walls and obstacles are black and the cells become lighter as long as the distance to these obstacles increases.

Implementation of a Smooth Slowness Potential
This implementation starts with the calculation of the logarithm of the inverse of the EVT of the 2D environment updated map (a priori + sensor data) (or the inverse of the EVT in the case of 3D maps). Each white point of the initial image (which represents free cells in the map) is associated to a level of grey that is the logarithm of the inverse of the 2D distance to the nearest obstacles. As a result of this process, a potential is obtained, proportional to the distance to the nearest obstacles, for each cell.
This function introduces a potential similar to a repulsive electric potential (in 2D), which can be expressed as � � � � ���� ���� � � � , where � represents the distance to the charge and � � , � � are constants. If � � �, the potential is � � These expressions of the potential � correspond to the electric potentials in 2D and higher dimensional cases.

Fundamentals of the Method
Maxwell's laws govern the propagation of electromagnetic waves. We can use the most typical simplification of the problem in isotropic and possibly non-homogeneous media for a monochromatic wave, considering the so-called almost flat waves. In this case, the optical wave propagates at wavelengths much smaller than image objects, and therefore the ray optics approximates wave optics. These equations, along with the mentioned approach, allow us to develop the theories based on rays, such as the geometrical optics, the theory of sound waves, etc. The Eikonal equation is derived from these equations. In the Sethian [17] notation ��������� � ������ (2) where ���� represents the distance to the initial set, ���� is the refraction index of the medium, and � � ��, �� in 2D or � � ��, �, �� in 3D. In geometrical optics, Fermat's least time principle for light propagation in a medium with space varying refractive index ���� is equivalent to the Eikonal equation. The Eikonal solution ���� is a scalar function whose isolevel contours are normal to the light rays, and the refractive index of a medium is the quotient of the light velocity in the vacuum and its velocity in the given medium. This equation is also known as the Fundamental Equation of the geometrical optics.
The most important aspect of this equation, from a path planning point of view, is that the solution is an exponential function of the refractive index where � � �� � is the angular frequency, which is constant because the considered wave is monochromatic, � is the period, ���� � � ���� is the refraction index, � is the vacuum light velocity, � is the wavelength, ���� is the anisotropic medium light velocity, and � � � � . In this expression, � and � are constants and the refraction index is a function of the position, i.e., the considered medium is anisotropic.
Since this solution is exponential, if the refraction index or first potential ���� is � � then the second potential ���� is also � � , and therefore the trajectories calculated by the gradient method over this potential willbe of the same class.

Robot Formations and Fast Marching
Robot formations motion planning deals with the problem of how to find the paths and postures between robots so that the whole formation can adapt to obstacles and other environment characteristics.
In this work, we have considered two types of formations: 1) leader-following, where two robots (followers) follow the leader as if they were bodyguards, forming a triangle among the three of them; 2) leaderprotecting, a formation of six robots in hexagonal configuration protecting the centre of the formation (leader). The leader-protecting leader will not be represented for the sake of clarity and simplicity.
As already stated, the FM method is based on a potential function without local minima that provides smooth trajectories. The main problem to deal with is how to apply this method to the robot formation motion, keeping its good properties. To solve this problem a repulsive force is needed between the formation robots so that they keep a security distance and do not end up crashing into each other. Furthermore, an internal attractive force is also required for maintaining the formation while adapting to the environment. The advantage of using the VFM method, as proposed here, is that each robot, at each time, has one single potential which is attracted to the objective and repulsed from walls, obstacles, and the other robots of the formation, while keeping the desired nominal inter-robot distances.
Our first approach to implementing the repulsive force was to add Gaussian functions to the distance potential of each robot, centred with the other formation robots, so as to keep the desired inter-robot distances balanced with the need to avoid obstacles and reach the goal. This approach is fast and works for a wide range of cases, though the behaviour is not natural in some special situations such as those presented in Figure 7 and Figure  8, for the corridor at the bottom of the images.
Taking this into account, the best solution to maintain distance between a robot and others is to set the robots as obstacles in the environment map, expressed in a matrix . The computational cost of this method is higher, since similar matrices , and (see next section) have to be computed at each step. However, this approach preserves the good properties of FM and avoids the appearance of local minima.
The solution adopted for the implementation of the repulsive force is detailed next and summarized in Figure 3.

Description of the Algorithm of Robot Formations Motion using Fast Marching
The objective of this algorithm is to integrate, in a potential without local minima, a force that attracts the robots to the goal, a force that repels the robots from the obstacles and a force that maintains the formation.

Base Method
As shown in the previous section, the VFM method is based on the construction of two potentials: the repulsive or first potential , which can be interpreted as the metrics, the viscosity, or the refraction index of the medium, and the second potential with a funnel shape, which represents the distance to the goal point measured with the metrics . The VFM algorithm consists of an initialization and a loop. In the initialization:


The initial map of the environment is read.  The path of the leader is calculated using VFM applied to this map.
Each cycle of the loop represents a step in the movement of the robots. Within the loop: 1. For each cycle , each follower robot has its own map of the environment, which is the initial map, with the other robots of the formation placed as obstacles: The map is a function of the global map and the local map , and the local map is the union (occupancy map) of the map sensed by robot sensors and the map with the other robots' positions: This matrix is different for each follower robot : This matrix is binary: its entries are 0 (black) in the wall and obstacle locations (including the other vehicles) and 1 (white) in the free space, but also it could be the probability of occupancy. Having the other vehicles as obstacles introduces a repulsive force between these vehicles. 2. The metrics matrix for robot is calculated applying the EVT operator to the matrix for each vehicle and then taking the logarithm of the result in order to have a repulsive field similar to the electric one. This matrix is obtained using the obstacles and walls initially set for the FM method. This way, the matrix is composed of real numbers, with value 0 in the walls and obstacles and real grey levels in the rest. The grey levels are darker near walls and obstacles.
The matrix represents the metric , the viscosity, or the inverse of the refraction in the Eikonal equation. 3. The partial goal points are calculated using the leader path and the desired formation. The distance from the partial goals to the leader path is proportional to the grey level of the partial goals position. This way, the formation tends to be near the initial positions and partial goals act as attractive forces between the vehicles. 4. The distance matrices are calculated applying the FM method to the metrics matrix for each robot, using as a goal the partial goal of the previous step.
5. The path of minimum distance measured with the metrics from each vehicle to its partial goal is found using its potential and the gradient method. 6. All the robots move a fixed number of points on the corresponding path. 7. The leader advances its path position a fixed number of points on the initialization path.
The aforementioned attractive and repulsive forces between the vehicles act as glue in order to maintain the formation, but allowing enough freedom to take the other obstacles, bends, and narrow passages into account. Some interesting aspects of the process are shown in Figure 5. Figure 5. Consecutive steps of the robots' formation travelling around the maze using the main method.

Adding Springs
In highly symmetrical situations, i.e., the robots have to traverse a small door and are initially placed orthogonally to the wall, the two follower robots have problems in passing together through the door. In order to solve this problem, it is necessary to introduce a precedence order. For example, if there is no room for the two robots to pass, the second one has to pass first, and then the third (the first is the leader).
The solution proposed to this problem is to calculate first the complete trajectory for the second follower and then the complete trajectory for the last follower. This provides an effect similar to using springs with different stiffness between the robots. This effect can be seen in Figures 7  and 9.
Another possible solution to solve these highly symmetrical situations is to find the VFM joined trajectory for two followers in the four-dimensional configuration space resulting from the union of the two followers.

Maximum Energy Configuration
Another interesting problem relates to what the formation must do if the nearest passage is narrower than the width of the robot formation, where the width of the formation is the orthogonal diameter to the movement direction. The solution is to specify the smallest possible formation size (or the maximum energy configuration) beforehand and to search for another possible passage larger than this size.
The solution proposed is to dilate the walls and obstacles with the minimum radius of the robot formation. Then, the trajectory of the leader is found using this dilated map. This way, it is ensured that the formation passes through the passage used by the leader trajectory. Using this trajectory, the rest of the algorithm is similar to the main method: 1. The partial goal points are calculated using the leader path and the desired formation. The distance from the partial goals to the leader path is proportional to the grey level of the partial goals' position. 2. The distance matrices are calculated using the metrics matrix for each robot, using as a goal the partial goal from the previous step. 3. The minimum distance path to the partial goal is found using the potential . 4. All the robots move a fixed number of steps on the corresponding path. 5. The leader advances its path position a fixed number of steps.

Reduction of the operational cost
Although the FM method is very fast (the computation time is 0.2 sec for the maps shown in Figures 6 and 7, which have a size of 628x420 cells), the proposed algorithm for robot formations has to calculate the FM wave propagation for each follower robot at each cycle of the algorithm. For this reason certain techniques should be used to make the algorithm faster. Since the FM method can be considered as the continuous version of Dijkstra's algorithm, our goal is to turn it into an almost one-dimensional algorithm. To achieve this, the VFM method is applied in a tube around the trajectory calculated for the leader. Thus, the propagation of the FM wave across the map is calculated only the first time to find the trajectory for the leader; the other times (once per cycle and follower robot) it is calculated in a tube around this trajectory, which drastically reduces the computation time. The steps are: 1. Enlargement of the trajectory calculated for the leader to obtain a tube. To achieve this, this trajectory is dilated to obtain a tube and the intersection between this tube and the map obtained from the environment (walls and obstacles) is used as the starting matrix (see Figure 9). 2. VFM-1st step. Using the map obtained in the previous step, a wave is propagated starting from the points representing the obstacles and walls. This is done with the Extended Voronoi Transform (also called Distance Transform in Computer Vision). The result is a potential map, which can be interpreted as a velocity map (or slowness map), as shown in Figure 9c. 3. VFM-2nd step. Based on the previous slowness map, the FM method creates a new potential that represents the propagation of the electromagnetic wave from the goal to the robot position. An even more important reduction of the computation time is related to the matrix t i D , as follows. As the vehicles are very close to the partial objectives, the expansion of the wave over the whole map is not necessary, but just its expansion into the tube, from the partial goal point to the corresponding vehicle. With this change the computation time of the matrix t i D , the most time-consuming part of the algorithm, goes from 0.2 sec to 0.016 sec, which implies an algorithm cycle of about 0.5 sec without parallelization (using a MacBook Pro platform at 3.06 GHz). The parallelized version of the algorithm has an algorithm cycle of about 0.3 sec and permits the use of many followers without increasing the computation time.

Simulation Results
The algorithms in the previous sections have been tested using Matlab simulations. A program has been created in which the initial and goal points are given and the VFM method is used for the leader and the rest of the formation. The method uses a nominal desired formation of the partial goals with a distance orthogonal to the leader trajectory proportional to the grey level of the actual position. The followers' trajectories to the partial goals are calculated using the VFM method. The results are shown in Figures 5 to 7.
The method has been designed for holonomic robots, but it is possible to apply these techniques for non-holonomic robots in [18].
The method has been used giving a sequence of points, joining these points with lines to obtain the leader trajectory. Using this trajectory, the followers' trajectories are calculated in the same way as in the main proposed method. In this case, the direction in the given points changes abruptly. It is important to study the behaviour of the method in these points, and we can observe that it behaves very well and even smoothes the followers' trajectories in comparison with the leader's, as depicted above.
The desired inter-vehicle distance has a maximum value in open areas and is proportional to the refraction index in the rest. This way, in small corridors the followers can be near each other. Figure 5 shows an example of a team of two vehicles following a moving leader using the main proposed method. The lines connecting each of the three vehicles represent the formation links between them. The circles represent the partial objectives that change at each step of the algorithm. The lines from the vehicles to the circles are the partial paths calculated using VFM. Figure 6 shows a similar case using the proposed method with springs of different stiffness. With this improvement the symmetry of the followers' paths is broken and it is easier for them to pass through small passages. Figure 7 shows consecutive steps of the robots' formation travelling around the maze using the main method, with a formation of six robots protecting the centre of the formation. The EVT can be considered as a viscosity or the inverse of the velocity of the medium. This way the method gives a maximum velocity in each point of the trajectory as shown in Figure 7e. Figure 8 shows consecutive steps of the robots' formation using the maximum energy configuration with springs. The formation does not use passages narrower than its maximum energy configuration.
In the simulations, the algorithms always determine a safe path where the robots avoid obstacles, keep the formation geometry in open space, and deform it when necessary to handle the presence of obstacles. Simulation results show that the method does not have local minima, which naturally results from the method used (light propagation has no local minima) and is a distinctive property with respect to potential fields and other similar methods in the class of motion planning algorithms. Furthermore, the method naturally ensures trajectory tracking by all the formation robots, since it naturally induces the required velocities at each point to reach the goal in minimum time, taking into account the constraints imposed by the presence of obstacles and the formation geometry.

Conclusions and Future Work
This paper presents a new methodology for the motion of a formation of holonomic robots. The formation is maintained by calculating the trajectory of the formation leader. In each iteration of the algorithm, using this trajectory the partial robot objectives are calculated. These partial objectives maintain the formation and have a variable distance between them proportional to the light velocity in that point (inverse of the refraction index). Each robot has an environment map with the other robots as objects. This map is used to construct the refraction index map or metric ( ) P x using the FM method. Using these metrics, the distance function ( ) D x representing the distance function measured with the metrics ( ) P x is built with the FM method. The partial trajectory of each robot is calculated with the gradient method.
The general trajectories have a behaviour like the light trajectory in a space with larger refraction index near the obstacles and walls and with an attraction force that tries to maintain the formation. The main contribution of the method is that it robustly achieves smooth and safe motion plans in real time that can be used at low control levels without an additional smoothing interpolation process. This allows the method to fuse collision avoidance and global planning in only one module, which can simplify the control architecture of the mobile robot, and without local minima, as in the case of the potential fields' original method and some other algorithms in the same class.
Additionally, it represents an innovation in formation control since it introduces an algorithm that simultaneously enables deformable formations, as in [8], and avoids the local minima problem of potential fields [14] by using a distance function based on a metric built with the FM method. This is particularly relevant when the formation needs to handle concave obstacles without some of the formation robots becoming trapped due to local minima of the inter-robot composition of attractive and repulsive potentials, as reported in [10]. Moreover, the path of the leader to the goal avoids obstacles and is free of local minima, since it relies on the VFM method.
Other local-minima-free methods could be used for the latter purpose, e.g. [19].
The future work relates to the introduction of sensor noise, uncertainty in the map (with connections with SLAM) and the introduction of small obstacles and moving obstacles.
All these questions can be implemented by adding Gaussian functions that model the uncertainty to the distance potential t i