Over the past decade, the bee population has been decreasing. One possible solution to this problem is the robotic bee. My project finds the optimal tour for a bee to navigate from a hive to a set of flowers and back to the hive again, occasionally also returning to the hive during the middle of the trip to recharge to avoid running out of battery due to its limitation of traveling at most distance m between charges. I solved this generalization of the traveling salesman problem with a function featuring both a slower-running method to find the shortest possible legal tour for the bee as well as a faster-running heuristic method that may return a suboptimal tour.

Introduction

Pollination is incredibly important in order for us to have a healthy thriving planet. More than 80 percent of the world’s flowering plants need a pollinator to reproduce, and bees do the lion’s share of that work. However, due to issues like habitat loss, pesticide use, and diseases, the bee population is decreasing. From April 2022 to April 2023, beekeepers across the United States lost 48.2% of their honey bee colonies. As the number of bees decrease, it could cause issues like food insecurity, lack of biodiversity, and potentially even extinctions. However, robotic bees could help solve this problem until real bees return. My project focuses on finding the optimal tour to pollinate a field of flowers for a robotic bee that can only travel a certain distance before it has to go back to its “hive” to charge.
First let’s find the coordinates of our flowers:
In[]:=
flowers=
;
Separate colors:
In[]:=
b=ColorSeparate[ImageResize[flowers,300]][[3]]
Out[]=
Binarize:
In[]:=
binarized=FillingTransform[DeleteSmallComponents[Binarize[b]]]
Out[]=
Find coordinates for ten flowers:
In[]:=
flowerLocations=Take[ImageKeypoints[binarized,KeypointStrength->0.015],10]
Out[]=
{{208.708,137.293},{182.193,89.7605},{98.703,36.9862},{162.096,123.074},{259.956,150.665},{276.807,112.19},{112.517,157.265},{252.727,125.495},{138.648,83.3644},{157.184,142.711}}

Finding Bee Tours

Let’s set the hive at the origin, and our maximum distance to 500:
In[]:=
hiveLocation={0,0};​​maxDistance=750;
If any flower farther from the hive than half the maximum distance allowed per trip from the hive, then no journey is possible:
In[]:=
If[AnyTrue[flowerLocations,EuclideanDistance[hiveLocation,#]>maxDistance/2&],Return[{}]]
Count the flowers:
In[]:=
flowerCount=Length@flowerLocations
Out[]=
10
Associate the index of each flower (starting at 2 since the hive is the first index) to its positional coordinates:
In[]:=
flowerPositionAssociation=Association@@MapIndexed[#2[[1]]+1->#1&,flowerLocations]
Out[]=
2{208.708,137.293},3{182.193,89.7605},4{98.703,36.9862},5{162.096,123.074},6{259.956,150.665},7{276.807,112.19},8{112.517,157.265},9{252.727,125.495},10{138.648,83.3644},11{157.184,142.711}

Finding the Shortest Hamiltonian Cycle

Each tour consists of a set of cycles each of which starts at the hive and visits a subset of the flowers before returning to the hive. Hence we desire a function to find the shortest Hamiltonian cycle through any complete graph whose vertices are the hive plus some subset of flowers.
Let’s look at the subset consisting of the first five flowers by index:
In[]:=
flowerIndices=Range[5]+1;
Create a distance matrix of the distances between the flowers, as well as between each flower and two copies of the hive (two hives in the same location are necessary because FindShortestTour[] requires starting and ending points to be different):
In[]:=
dm=DistanceMatrix[Join[flowerPositionAssociation/@flowerIndices,{hiveLocation,hiveLocation}]]/.0.->Infinity
Out[]=
{{∞,54.4273,148.87,48.7322,52.964,249.816,249.816},{54.4273,∞,98.7713,38.906,98.7742,203.104,203.104},{148.87,98.7713,∞,106.91,197.295,105.405,105.405},{48.7322,38.906,106.91,∞,101.675,203.524,203.524},{52.964,98.7742,197.295,101.675,∞,300.461,300.461},{249.816,203.104,105.405,203.524,300.461,∞,∞},{249.816,203.104,105.405,203.524,300.461,∞,∞}}
Create a weighted graph from distance matrix:
In[]:=
gra=WeightedAdjacencyGraph[dm,VertexLabels->"Name",EdgeLabels->"EdgeWeight"]
Out[]=
Find the shortest tour that starts at one hive, visits all the flowers in the subset, and returns to the duplicate hive:
In[]:=
{dist,tour}=FindShortestTour[gra,VertexList[gra][[-2]],VertexList[gra][[-1]]]
Out[]=
{608.171,{6,3,2,5,1,4,7}}
Drop the hives from the beginning and end of the tour, at return the flowers to their original indices:
In[]:=
flowerIndices[[tour[[2;;-2]]]]
Out[]=
{4,3,6,2,5}
Here is the full code for the function:
In[]:=
findShortestHamiltonianCycle[hive_,flowerPositionAssociation_][flowerIndices_]:=​​Module[{dm,gra,dist,tour},​​dm=DistanceMatrix[Join[flowerPositionAssociation/@flowerIndices,{hive,hive}]]/.0.->Infinity;​​gra=WeightedAdjacencyGraph[dm];​​{dist,tour}=FindShortestTour[gra,VertexList[gra][[-2]],VertexList[gra][[-1]]];​​{dist,flowerIndices[[tour[[2;;-2]]]]}​​]

Optimal Bee Tour

The optimal tour can be found by considering all partitions of the set of flowers, then finding the one that admits the shortest legal bee tour with stops to the hive between tours of each partition part.
Associate subsets of flower indices to lists whose first element is the total distance required to travel from the hive, visit all the flowers in the subset, and return to the hive (ignoring the maximum distance constraint) and whose second element is the list of indices of flowers in the subset ordered by how they are visited in the shortest journey:
In[]:=
subsetDistanceAssociation=AssociationMap[findShortestHamiltonianCycle[hiveLocation,flowerPositionAssociation],DeleteCases[Subsets[Range[Length@flowerLocations]+1],{}]]
Out[]=
{2}{499.632,{2}},{3}{406.209,{3}},{4}{210.81,{4}},{5}{407.049,{5}},
⋯1016⋯
,{2,4,5,6,7,8,9,10,11}{718.028,{8,11,2,6,7,9,5,10,4}},{3,4,5,6,7,8,9,10,11}{721.477,{4,10,3,9,7,6,5,11,8}},{2,3,4,5,6,7,8,9,10,11}{721.498,{8,11,5,2,6,7,9,3,10,4}}
Full expression not available
(
original memory size:
498.2 kB)
Find all partitions of the set of flower indices such that no part in the partition violates the maximum distance condition:
In[]:=
Quiet[<<Combinatorica`]​​flowerPartitions=Select[Combinatorica`SetPartitions[Range[flowerCount]+1],Max[(subsetDistanceAssociation/@#)[[All,1]]]<=maxDistance&]
Out[]=
{{{2,3,4,5,6,7,8,9,10,11}},{{2},{3,4,5,6,7,8,9,10,11}},{{2,3},{4,5,6,7,8,9,10,11}},{{2,4,5,6,7,8,9,10,11},{3}},{{2,3,4},{5,6,7,8,9,10,11}},
⋯115966⋯
,{{2,9},{3},{4},{5},{6},{7},{8},{10},{11}},{{2,10},{3},{4},{5},{6},{7},{8},{9},{11}},{{2,11},{3},{4},{5},{6},{7},{8},{9},{10}},{{2},{3},{4},{5},{6},{7},{8},{9},{10},{11}}}
Full expression not available
(
original memory size:
60.4 MB)
Find the partition that will give the shortest tour distance:
In[]:=
flowerPartitions=First@MinimalBy[flowerPartitions,Total[(subsetDistanceAssociation/@#)[[All,1]]]&]
Out[]=
{{2,3,4,5,6,7,8,9,10,11}}
Replace each part of the flower partition with the list of its distance and its optimal ordering:
In[]:=
flowerPartitions=subsetDistanceAssociation/@flowerPartitions
Out[]=
{{721.498,{8,11,5,2,6,7,9,3,10,4}}}
Total the distances and combine the flower tours, riffling in the hive’s index 1 between each journey through flowers:
In[]:=
flowerPartitions={Total[flowerPartitions[[All,1]]],Join@@Riffle[flowerPartitions[[All,2]],{{1}},{1,-1,2}]}
Out[]=
{721.498,{1,8,11,5,2,6,7,9,3,10,4,1}}

Heuristic Method

The number of nonempty subsets of a set of
n
flowers is
n
2
-1
, and the number of partitions of the complete set is given by the
n
th Bell Number, which we see can grow quite large:
In[]:=
Table[BellB[n],{n,20}]
Out[]=
{1,2,5,15,52,203,877,4140,21147,115975,678570,4213597,27644437,190899322,1382958545,10480142147,82864869804,682076806159,5832742205057,51724158235372}
Hence for large flower sets we desire a method faster than the one detailed above, even if it may fail to find the optimal tour for the bee. This heuristic method will use the function FindClusters[] to cluster the flowers:
Find the shortest tour of each of these flower clusters:
If the tour does not violate the maximum distance condition, return it; otherwise, keep incrementing the number of clusters until a legal tour is found, and return that. Along the way, store the shortest Hamiltonian cycles of each cluster so that it won’t have to be computed again. Details of implementation can be found in the full code below.

Full Code

Using the Function

Set the hive at {1,1}, pick ten random points for the flowers, and set the max distance at 5:
Find the optimal tour:
Find a tour using the heuristic method:
Watch the bee make its rounds:

Conclusion

This project investigated how to find the shortest tour a robotic bee could take when pollinating a field of flowers, given the constraint that the robotic bee would have to go back to the “hive” to charge before reaching a certain maximum distance since the previous charge. We provided a function capable of producing an optimal tour through brute force as well as a potentially suboptimal tour using a heuristic method. In the future, I’d like to improve the correctness of the heuristic model while still maintaining the speed.

References

Randall, B. (2020, June 22). The Value of Birds and Bees. Farmers.gov; USDA. https://www.farmers.gov/blog/value-birds-and-bees​
Planet Bee Foundation. (n.d.). Why Bees are Dying | Resources | Planet Bee Foundation. Planet Bee. https://www.planetbee.org/resources/why-bees-are-dying

Acknowledgements

Thank you to my mentor Dan McDonald for supporting me on my project. This project would have been significantly more difficult to complete without his help.

CITE THIS NOTEBOOK

Tour optimization for robotic bee pollination​
by Grace Anghelescu​
Wolfram Community, STAFF PICKS, July 10, 2025
​https://community.wolfram.com/groups/-/m/t/3498845