MULTI-AGENT PATH PLANNING OF ROBOTIC SWARMS IN AGRICULTURAL FIELDS

: Collaborative swarms of robots/UAVs constitute a promising solution for precision agriculture and for automatizing agricultural processes. Since agricultural ﬁelds have complex topologies and different constraints, the problem of optimized path routing of these swarms is important to be tackled. Hence, this paper deals with the problem of optimizing path routing for a swarm of ground robots and UAVs in different popular topologies of agricultural ﬁelds. Four algorithms (Nearest Neighbour based on K-means clustering, Christoﬁdes, Ant Colony Optimisation and Bellman-Held-Karp) are applied on various farm types commonly found around Europe. The results indicate that the problem of path planning and the corresponding algorithm to use, are sensitive to the ﬁeld topology and to the number of agents in the swarm.


INTRODUCTION
By 2050, nearly 10 billion people will live on the planet. Looking at the population predictions, the first question which comes to mind is: "is there enough food for our future?" Unfortunately, current worldwide crop yields (even with its existing increasing trends) would not satisfy this increasing food need. Based on the current crop yields, researchers suggest production of agricultural yields must be increased by roughly 50% (Article, 2019). In order to close this production gap, using high-tech and automated precision farming systems, as well as approaches borrowed from organic farming have seen as the ultimate way which could boost yields in the existing farming places (Article, 2019). Therefore, researchers have started to focus on finding the best ways to use intelligent systems and robotics for overcoming the food security problem in times of climate change (Balafoutis et al., 2017). Precision agriculture is nowadays the name given to technological developments for assessing and responding to local variations and emergencies in crop health on farms. It involves collecting, interpreting and acting on data which has potential to be automatized with robotic systems. A popular metric to learn the vegetation health is to look at the normalized difference vegetation index (NDVI) (Weier and Herring, 2000). It is calculated from images and indicates whether a region contains healthy vegetation or not. Satellites or robots (especially unmanned aerial vehicles, called UAVs (Hassler and Baysal-Gurel, 2019)) are used to collect such images over farm lands. A current research area is improving this automation, in particular managing robot swarms for precision agriculture (Carbone et al., 2018). This includes automatic design of ways that UAVs do exploration, collaboration on tasks, communicate and find the optimal fleet formations. Since ground robots have limitation with mobility (not being able to go over high vegetation and collect images of the area from a larger viewing point), the most appealing automation solution for agricultural fields have been based on UAVs. Therefore, some remote sensing researchers have already proposed algorithms for UAV based crop monitoring (Valente et al., 2019), (Finn et al., 2019), * For communication, email address of the corresponding author: a.kamilaris@rise.org.cy (D'Urso et al., 2018), (Honkavaara et al., 2012), (Ok and Ok, 2017). However, usage of multiple UAVs (UAV swarms) is discussed in very limited number of studies (Goel et al., 2016), especially in relation to agricultural applications (Hassler and Baysal-Gurel, 2019). This is mainly because of difficulty of path planning for multiple robots who are trying to find the optimal paths for themselves in order to work collaboratively on a field. Since this multiple-agent path planning field is still open for intelligent algorithms to guide UAVs on different farm lands, herein we focus on an aspect of collective exploration, path planning through agriculture fields. That is, allocating crops and assigning paths for each UAV (robot) to follow, in order to complete the tasks in an efficient way as a team. Two main scenarios are considered. In the first, the goal is to plan paths through the entire field. In literature, this is known as the coverage path planning (CPP) problem (Oksanen and Visala, 2009b). The second scenario is to plan paths through only selected points or areas. These points could be preassigned by the farmer or ground workers, or from an image from a satellite or a single high altitude UAV. This scenario is closely related to extensive work on the vehicle routing problem (VRP) in operations research (Laporte, 1992a), (Laporte, 1992b), (K. Braekers and Nieuwenhuyse, 2016). These path planning problems are NP-hard. This means that the solution space grows very rapidly with the problem size, and there is no efficient algorithm for finding an exact optimal solution. Therefore in general it is best to search for an approximate best solution. The focus then was on investigating algorithms for path planning. Two models were implemented: a continuous world model (no obstacles, for UAVs), and a grid world with obstacles (with obstacles, for ground robots).

RELATED WORK
Swarm robotics is the study of controlling multiple autonomous robots separately but achieving a desired collective behaviour (Brambilla et al., 2013). These behaviours include flocking together while travelling, working together to complete a task, or collective exploration. The particular task considered in this paper, path ISPRS Annals of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume V-1-2020, 2020XXIV ISPRS Congress (2020 planning, is part of collective exploration. Swarm robotic is a relatively new possibility for agriculture with much potential (Carbone et al., 2018). The aim of path planning is to explore an entire area as efficiently as possible (Oksanen and Visala, 2009a). It usually has the additional constraints of simple trajectories to simplify control for autonomous robots. It is a NP-hard problem. Moreover, there are additional challenges when it is applied to ground vehicles in agriculture. For example, a common constraint is that large vehicles can only turn in the headland areas (Hameed et al., 2010), i.e. these are areas on the edges of the field, outside of the main crop regions. Generally, for covering an entire farm, there is an added complexity due to the obstacles inside the farm and the irregular shapes of farms (Machl et al., 2013). We tried to address some of the aspects during problem modelling (see Section 3.1).
Regarding path planning of multiple field robots, Seyyedhasani and Dvorak (Seyyedhasani and Dvorak, 2017) converted the path planning problem to a VRP. Then, they used a modified Clarke-Wright algorithm to find a non-optimized solution, and 600-700 iterations of a tabu search algorithm to optimise this solution.
Regarding path planning of multiple UAVs in agricultural fields, there is very limited research, perhaps due to the fact that UAVs have relatively recently entered the agricultural area. Barrientos et al. (Barrientos et al., 2011) implemented a path planning for a small fleet of 3 UAVs in a 6.4ha farm. Their path planning algorithm transformed the farm to a low-resolution two-dimensional grid, divided the area equally in between the agents, and then used a search algorithm with backtracking to find paths that minimised the number of turns. Ju and Son (Ju and Son, 2018) applied a distributed swarm control algorithm to suggest that the performance of the multi-UAV system is significantly superior to the single-UAV system, considering seven performance metrics: total time, setup time, flight time, battery consumption, inaccuracy of land, haptic control effort, and coverage ratio. Albani et al. (Albani et al., 2017b), (Albani et al., 2017a) investigated a multi-agent system for field coverage and weed mapping, based on a reinforced random walk with inhibition of return, where the information available from other agents is exploited to bias the individual motion pattern. A problem with the aforementioned papers (and the algorithms used by the authors) is that they provide very constrained solutions. The generated paths cannot navigate around complex geometries or they require the fields to be rectangular or made out of simple polygons (Seyyedhasani and Dvorak, 2017), (Ju and Son, 2018), (Albani et al., 2017b), (Albani et al., 2017a). Some of the single-agent solutions cannot generalize to the multiple-agents case (Bochtis and Zhou, 2015), (Oksanen and Visala, 2009a), (Hameed et al., 2010) or the solutions have been tested on a few agents only (Barrientos et al., 2011), (Ju and Son, 2018) (i.e up to 3 UAVs). Finally, all aforementioned related work considered one or two algorithms, but did not try to use and compare a variety of possible solutions to the problems under study.

METHODOLOGY
From this point on, the word agent will represent both ground robots and UAVs in the modelling of the problem. These words will be used interchangeably.

Problem modelling
The goal is to plan paths for multiple agents through interest points spread over a large area. These agents can either be aerial or ground vehicles. Two two-dimensional (2D) environments were created: a polygon-shape structure for the aerial vehicles and a grid structure for the ground vehicles.

Modelling of agricultural fields
The agricultural field structures for the ground robots and for the UAVs are depicted in Figures 1 and 2 respectively. For the ground robots, the fields were designed considering the most popular structures based on equally popular crops. For example, Figure 1 (left) is representative of olive trees, citrus fruit trees, almond trees and generally fruit and nut trees. Figure 1 (middle) is representative of various vegetables such as pepper, tomatoes, lettuce and also of legumes, cornfields, cereals, soybeans, cornflower etc. Figure 1 (right) constitutes another typical example of a field of fruit trees, in which some obstacles exist (e.g. house, irrigation system, greenhouse, storage room etc.). For the UAV case, a different design logic has been considered. UAVs usually need to fly above specific regions of the fields where some abnormality occurs (e.g. water stress, lack of enough fertilizer, pests have been detected, areas with low NDVI index etc.). A likely case of such abnormality in some field is illustrated in Figure 2 (left). To capture this function of UAVs, which is currently the most typical and popular one in agri-fields, we used points randomly distributed in small clusters (see Figure 2, middle to right). We refer to these points as points of interest (POIs). There was no need for placing the crops in regular patterns because on one hand the selected regions would be irregular in real-life situations and, on the other hand, the UAVs are not constrained by paths on the ground. The number of POIs was chosen to be 300 for each field for the UAV case.
In the case of the ground robots, obstacles were implemented, each occupying one or more cells (small and continuous lines in Figure 1). Obstacles represent crops or areas of the field where the agent is not allowed to pass. Obstacles are not necessary in the UAV case, as they fly above the crops. Simulations run until all POIs are visited, while return journeys back to the basis are estimated afterwards. In particular, for the fields in Figure 1, we assume a number of POIs varying from 100 to 150. Each point at the 2D plane has a unique position (x, y). Agents are modelled as point particles too. Their state is given by their current 2D position (x, y) and their linear velocity υ, which is uniform for all agents. The agents always started at position (0, 0), where they had to return afterwards. Each iteration of the algorithms implemented triggers an increment of a time step ∆t. The next agent state is calculated as: 3.1.2 Cost function The cost function used as the optimization goal of the path planning algorithms (see Section 3.3 is the one suggested by Seyyedhasani and Dvorak (Seyyedhasani and Dvorak, 2017) and it is given below: where d k is the distance travelled by agent k and m is the total number of agents. Similar to (Seyyedhasani and Dvorak, 2017), T is the maximum distance while D is the aggregated distance travelled by all agents. Thus, the cost function is defined as the weighted combination of maximum distance travelled and average distance travelled by each agents. Parameter z defines the weight of maximum distance vs. average distance, allowing to balance these two requirements towards an optimal solution that considers both total efficiency and fairness of individual agents.
ISPRS Annals of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume V-1-2020, 2020 XXIV ISPRS Congress (2020 edition)  Cost function has been used mainly for the ant colony optimisation algorithm (see Section 3.3).

Approximation to real fields
In order to approximate real-world agricultural field areas and distances, the following decisions were made: • For the UAV case, the velocity υ was set to 10m/s or 36km/h. This is a realistic speed for real UAVs (Carbone et al., 2018).
• For the ground robots case, velocity υ was set to 0, 83m/s or 4km/h, which is a realistic speed for ground robots.
Each iteration triggers an increment of a time step of 0.2s. In the context of this paper, velocities are used to estimate the transition of a robot/UAV from one grid cell to the next, according to one time step increment at each algorithm's iteration.
In theory, maximum distance T of the cost function can be used as a proxy for time, i.e. the agent that travels the longest distance will be the slowest agent. By knowing the average speed of the robot/UAV, one can easily then calculate actual times. This is an important aspect of future work (see Section 5.1). Each grid cell on the 2D plain represents an area of 2.25m 2 . The total area of each field shown in Figure 1 was around 0.3ha (hectares). The total area of each field of Figure 2 was around 20ha. These values match real-world agricultural farms.

Assumptions
The following assumptions were made to further confine the scope of the problem: 1. Number of agents: One up to five agents can be employed, either ground robots or UAVs. 2. Agent types: There is only one type of agent. 3. Origin: The agents start at a common location. 4. Destination: All agents must return to their start locations. 5. Scenario: Static, i.e. the field layout and POIs do not change, which is typical for an agricultural field on a typical day or even week. 6. Observability: Global, i.e. the field layout and all the POIs are known in advance. 7. Space: Agents can move freely in the two-dimensional space. 8. Interactions: The time needed for performing agriculturerelated tasks at each interest point (e.g. pruning, spraying, crop collection, taking photos etc.) is ignored (Hassler and Baysal-Gurel, 2019). Interactions between agents (i.e. collisions) are also ignored.
NN and Christofides constitute approximation algorithms. They are usually faster, as they have a limited number of iterations based on heuristics. In particular, NN was chosen for its simplicity and easy implementation. Christofides provides a bounded worst case solution and it is a very popular benchmark for solving VRP problems. However, they are not very flexible in their use generally. ACO is a meta-heuristic algorithm, i.e. it is by design more flexible (e.g. it can automatically choose the number of agents to employ), but also more unstable. It therefore requires careful tuning to produce good results. It is generally slower than 1 BHK was used as a comparison baseline for the run time of the other algorithms (Section 4.2) and it was not tested further on the different agricultural fields (Section 4.3 and 4.4).  NN and Christofides, while there is always the risk to get stuck in local minimum (missing the global, optimal minimum of the optimization function). BHK belongs to the family of dynamic programming algorithms. It is generally a highly accurate algorithm that is used mostly on small problem sets because of computational limitations. BHK has extremely high running complexity (O(n 2 2 n )) and space complexity (O( √ n2 n )). A comparison of strengths and weaknesses of the chosen algorithm is shown in Table 1.

RESULTS
In this section the results are presented. First, the choice of the parameters used in the ACO algorithm is discussed in Section 4.1, since ACO required much fine-tuning in order to work in the most optimal way. Then, the run time of the different algorithms, introduced in Section 3.3, are presented together with the choice of the candidate algorithms for the tests on the fields presented in Figure 1 and 2. Finally, the results obtained by applying the three proposed algorithms are shown, in Section 3.3. The algorithms are tested on three different fields using UAVs (see Figure 2) in Section 4.3 as well as on three different fields using ground robots (see Figure 1) in Section 4.4. As mentioned in Section 3.2, the maximum number of agents has been defined as 5, since it is assumed that 5 ground robots or UAVs are enough to effectively cover a complete field. The three algorithms are compared in terms of the following aspects: Among the candidate algorithms, mainly 2 ACO uses the cost function of Equation 2, thus we selected the metrics described above for fair comparison of the four algorithms. In Figures 3, 5 and 6a-11a, total distance is illustrated with a solid line, average distance with a dash line, while minimum/maximum distance travelled by agents is depicted by the shaded areas. Examples of actual solutions are shown in Figures 6b-11b. Table 2 summarizes the choice of values for the parameters information elicitation coefficient α, the expected heuristic coefficient β and the pheromone evaporation coefficient ρ, since this combination gave the best results after an exhaustive run of the algorithm combining all possible combination of values.

Parameter
Value Information elicitation coefficient α 3 Expected heuristic coefficient β 5 Pheromone evaporation coefficient ρ 0.01 Weight of cost function z 0.5 Table 2: Choice of parameters for ACO.
The second important decision refers to the parameter z in the cost function (see Section 3.1.2), which allows a trade-off between the two objectives of the equation. In Figure 3, the total distance (solid line), the average distance (dash line) and the min/max distance (shaded area) are shown with respect to the number of agents while varying the parameter z of the cost function. For z = 0, the first term of the cost function disappears, thus the only objective of ACO is to minimize the average distance (see Figure 3, blue line). This leads to one or more agents being idle, as can be seen in the figure (blue shaded area). When z = 1 (see Figure 3, orange line), the workload of each agent is well balanced leading though to a larger total distance D travelled. Based on these considerations, z = 0.5 (see Figure 3, green line) has been selected for all the following experiments because it offers a good trade-off between the minimization of the total distance D and the fair spreading of the workload among the agents.

Algorithms' run time
Run time of each algorithm is a crucial aspect because it defines the time of planning required for each field. This planning can change almost daily in modern precision agriculture scenarios. We measured the time needed to find an optimal solution with respect to the number of POIs (see Section 3.1). These tests were performed using a Intel i7-8550 CPU with 1.80 GHz and 8GB of RAM. In particular, we compared the following algorithms on finding the shortest path to visit all the POIs using UAVs: BHK computes the global optimal solution, but at the price of high computational/time and space complexity (see Section 3.3). As shown in Figure 4, the optimal solution computed with the BHK algorithm requires exponential time with respect to the number of POIs and agents. With approximately 20 POIs, the BHKbased solutions require more that 8GB of RAM and their run time is 10 times higher than ACO. Since the number of POIs was chosen to be 300 for UAVs and 100-150 for ground robots (see Section 3.), BHK was not used for the experiments presented in Sections 4.3 and 4.4, due to these expensive computational requirements.
Christofides and NN ranked algorithms (see Figure 4) require the shortest amount of time, while time does not increase exponentially with the number of points to visit. ACO requires more time than Christofides and NN ranked algorithms due to its higher complexity, but it follows a similar trend as the number of points grows. Eventually, the performance of the algorithms with the metrics described in Section 4. is shown in Figure 5. Even though NN has the shortest run time, it achieves the worst performance in the simple test case, shown in Figure 5 (right). Because of this reason, NN is excluded for the real experiments presented in Section 4.3. The candidate algorithms that will be presented in the following section are: Christofides' algorithm with the different clustering techniques presented in Section 3.3 and ACO with the parameters choice described in Section 4.1

Results for the UAVs case
For UAVs, we compared ACO, using the cost function defined in Section 3.1.2 and parameters as shown in Table 2, with two different clustering techniques for the Christofides algorithm: sector or rectangular clustering and K-mean clustering with Cartesian coordinates.
When the field topology is symmetric without constraints (i.e. obstacles) (see Figures 6b and 8b), Christofides algorithm employing both clustering techniques outperforms ACO in terms of minimization of the total and average distance travelled and also in terms of fair sharing of the workload among the agents (see Figures 6a and 8a respectively). When symmetry in the field topology is lacking (see Figure 7b), ACO achieves performance comparable to the Chistofides algorithm with K-means clustering, while it performs better than the Chistofides algorithm with sector/rectangular clustering (see Figure 7a).

Results for the ground robots case
For the ground robots case, we compared once more ACO, using the cost function defined in Section 3.1.2 and parameters shown in Table 2, with different clustering techniques for Christofides algorithm: sector or rectangular clustering and K-mean clustering with Cartesian and polar coordinates. In particular, we assumed squared obstacles as in the fields shown in Figure 1 (left and right), while we assumed line obstacles as in the field in Figure 1 (middle).
Similarly to the UAVs case, when the field presents symmetry (see Figures 9 and 11), K-means clustering with polar coordinates outperforms ACO. However, in contrast to the results of the UAVs case, as soon as the number of agents grows, ACO tends to improve its performance (see Figures 9a and 11a). In the field shown in Figure 1 (right), it achieves the best performance among the compared algorithms (see Figure 11). When the field is not symmetric and heavily constrained, as in Figure 10, ACO outperforms Christofides algorithm with both clustering techniques.

DISCUSSION
Our results indicate that different algorithms work differently and are suitable for specific agricultural scenarios of different topology and also of different type of agents (i.e. ground robots vs. UAVs). Table 3 summarizes our findings on the specific modelling choices we made. No obstacles refers to the setting for UAVs (see Figure 2), vertical lines and horizontal lines to fields with vertical/horizontal alignment of crops respectively (see Figure 1,middle) while squared setting is about fruit/nuts trees (see Figure 1, left). Vertical/horizontal lines are always defined in respect to the agents' initial placement. From Table 3, it seems that ACO works better when minimizing the total distance needed to be travelled by the agents (i.e. D in Equation 2) as well when minimizing the maximum distance that needs to be travelled, which corresponds to total time needed (i.e. T in Equation 2). This holds both for UAVs and ground robots. However, when a balanced solution is sought between distance/time (i.e. z = 0.5 in Equation 2), Christofides algorithm comes into play as a more efficient one, either together with Kmeans for squared and horizontal fields as well as for UAVs without boundaries in the fields, or together with sectors for ground robots and vertical lines at the fields. These results show that agricultural fields are sensitive to their topology and type of robot swarms used, in terms of the best algorithm to be used for path planning. It is also the case that the number of available agents affects the best algorithm to be used (see Figures 9a and 11a).
Although the scenarios considered are static in the sense that the agents know in advance the total plan, this can become dynamic (especially) in the case of UAVs whose POIs depend on daily observations made by remote sensing, towards precision agriculture. Thus, the optimal use of path planning algorithms becomes quite important and their running time affects the practicality of each solution (see Section 4.2). This is the reason BHK does not seem to be a good choice in such scenarios.   Table 3: Summary of findings.

Future work
As a future work, we will spend more time on implementation of other multi-agent path planning methods in order to make more fair judgement of our methods. Besides, we will focus on three mathematical issues of our own method: better optimization of the parameters (e.g. choice of z), quantitative analysis of the sensitivity of the algorithms to the fields' topology, testing them on various scenarios; the generalization of our method for different farm shapes automation of the progress (from imaging to path planning and robotic interaction with the crops), we will extract the target positions from NDVI maps with computer vision algorithms . We would like to develop an intelligent solution that, provided with the list of POIs to visit and the knowledge of field topology, automatically chooses the optimal number of agents to deploy and optimally plans their paths. Moreover, we will consider actual constraints of UAVs (e.g. max. flight time) and ground robots (e.g. max. distance travelled in a number of hours), adapting the cost function and z parameter in a way to reflect also time limitations. This will require the use of more realist physical models for UAVs and ground robots. Finally, a dynamic scenario will be explored, in which the agents do not have prior knowledge of the field but need to combine their ad hoc knowledge in order to address the optimization problem in the best way possible.

CONCLUSIONS
This paper has tackled the problem of optimizing path routing for a swarm of ground robots and UAVs in different topologies of agricultural fields. Four algorithms (Nearest Neighbour based on K-means clustering, Christofides, Ant Colony Optimisation and Bellman-Held-Karp) have been applied on various farm types which are commonly found around Europe. The results indicate that the problem of path planning and the corresponding algorithm to use, are sensitive to the field topology and to the number of agents in the swarm. We believe that the proposed method used in this paper provides important information towards developing multi-agent robotics, in order to optimize the time, distance and energy spent in the farms by the swarms of agents, towards more precision farming. Further experimentation in this open research field is still required.