Development of Heuristic Approaches for Last-Mile Delivery TSP with a Truck and Multiple Drones

.


Introduction
In recent years, last-mile delivery has witnessed an enormous increase in quality demand in a continuously growing market, especially due to the rapid spread of online shopping platforms.Usually, e-commerce retailers are required to deliver a huge number of parcels within small time windows, to several customers and possibly without charging any delivery fees.These requirements are obviously conflicting and they shall be tackled with efficient management.
By definition, last-mile delivery refers to the final stage of the logistical process, when the package is transported from the transportation hub to the destination of the final client.It follows the main stages of the last-mile delivery process: (1) when a consumer places an order, it is entered into a digital system where it can be tracked; (2) at the transportation hub, where orders are waiting to be delivered, the logistic process starts; (3) orders are strategically assigned to delivery employees based on delivery addresses with the aim of optimizing the process itself; (4) orders are scanned before being loaded onto delivery trucks in order to monitor the parcel on both the customer and corporate sides; (5) orders are delivered successfully to their intended location, and the delivery status is updated after a proof of accurate shipment is requested.
The fundamental paradox is also that, despite the fact that customers expect quick and possibly free delivery services, the last-mile delivery process is the most expensive one Comprehensive literature reviews of the applications of drones in the last-mile delivery domain are proposed in [14][15][16].An example of an application is represented by the multiple Flying Sidekick Travelling Salesman Problem (mFSTSP), where a fleet of UAVs is used in combination with a truck that simultaneously performs other delivery tasks and can also serve as a platform for the launch and the recovery of the drones.A similar application is considered in this work.
The problem of scheduling parcel delivery tasks to a truck operating in combination with a fleet of UAVs is a variant of the Travelling Salesman Problem (TSP).In computer science, the TSP is still an open problem to be solved in a globally optimal way due to the intractability of its inherent combinatorial nature.A problem is said to be intractable if the computational complexity is superpolynomial.The TSP is defined as an NP (Nondeterministic Polynomial) problem, which means that a solution can be verified in polynomial time, but the optimal solution cannot be computed in polynomial time.However, it has not been proven yet that a polynomial resolution algorithm does not exist for these kinds of problems.In this case, scheduling parcel delivery tasks to a delivery system made by a truck and many drones belongs to the class of NP-hard problems, i.e., problems hard at least as the hardest problem in the NP-class of decision problems.The fact that the considered TSP problem can be formulated as a decision problem is straightforward.Due to the intractability issue of the optimal solving algorithms of the mFSTSP, and any TSP problem in general, several optimization methods have been developed over the last decades to address NP-hard problems.They are divided mainly into four categories based on the process of inspiration: evolution (e.g., genetic algorithms [17]), hunting (e.g., particle swarm optimization [18,19]), nature science (e.g., simulated annealing [20]) and behavioral/strategic processes (e.g., colony search optimization [21,22]).Refer to [23] for a comprehensive survey about the TSP problem's applications, the state-of-the-art solution approaches, and the complete taxonomy.
In this paper, we propose a mathematical programming statement of a particular extension of the TSP: the mFSTSP.Trivially, this problem becomes increasingly difficult if realistic constraints are considered and added to the formulation, such as the modeling of the tasks of launching and retrieving the drones, the parcel service time to complete the delivery, the drones' battery discharge, the synchronization between the truck and drones, the presence of more than one drone in the fleet, etc.Thus, we believe our work provides a significant contribution to the current state-of-the-art about the mFSTSP formulations and solution methods, since our formulation accounts also for some last-mile delivery specific issues, such as the service delivery time of the parcels, which is an inefficiency factor of the process, as highlighted in [10].
Due to the NP-hard nature of the addressed problem, heuristics approaches can often provide efficient and effective solutions.Heuristics algorithms do not guarantee an optimal solution but can provide good approximations in a reasonable amount of time.On the contrary, commercial solvers may not be suitable to solve our mFSTSP problem due to the addition of constraints and variables associated with the drones and their interactions with the truck.The problem becomes more complex, especially as the number of drones and the size of the problem increases.As also discussed in [24,25], the computational requirements of finding an exact solution to the mFSTSP using commercial solvers can become impractical or time-consuming, even for small-scale instances.For these reasons, we implement and evaluate a few sub-optimal heuristic resolution approaches: a local search solution, an evolutionary-based solution with two variants, and a lightweight ad-hoc greedy solution.
The main contributions of this work are (i) the formulation of a mathematical programming representation of a well-defined delivery scheduling problem with a truck and multiple drones and considering numerous constraints, which are generally partially evaluated at the state-of-the-art, (ii) the design of a hybrid evolutionary-based approach to be compared to both a computationally fast greedy solution method and the well-known local search method, (iii) the presentation of realistic simulation results in an urban environment.This paper is organized as follows.Section 2 discusses the related work conducted in recent years concerning routing and task scheduling algorithms for last-mile delivery with truck and drones.Section 3 presents the problem formulation and the assumptions, as well as the proposed solution approaches and the simulation settings.Section 4 presents and discusses the simulation results obtained by comparing the performances of the proposed methods in terms of the total completion time of the delivery cycle for two scenarios of reference.Our conclusions and future research directions are drawn in Section 5.

Related Work
A literature review of heuristic routing solutions and mathematical formulations for tandem truck-drone problems is presented in the following.
The work in [26] proposes a mathematical formulation of TSP for an energy-efficient delivery system based on the use of UAVs in cooperation with the public transport system, showing the potential to optimize drones' energy consumption and battery charge by exploiting public transport vehicles as carriers.In [27], taking into account truck and drone's synchronization and minimization of total delivery time, two heuristic algorithms are proposed for solving large-size problems: Drone Truck Construction (DTC) and Large Neighborhood Search (LNS).The DTC algorithm tackles the routing problem of the first truck route from the depot to the scheduled customers, and the LNS algorithm schedules the drone routes considering the truck as a movable depot.The research in [28] suggests a hybrid genetic search algorithm based on a split algorithm, crossover, and local search operations.A restore strategy is designed to enhance the convergence properties of the algorithm, as well as an adaptive penalization mechanism to dynamically balance the search between feasible and infeasible solutions in order to solve the TSP-D (Travelling Salesman Problem with Drones), in the single drone variant.Simulation results show the capability of the evolutionary-based method to efficiently handle both min-cost and min-time optimization.Another hybrid genetic algorithm is proposed in [24] to tackle the mFSTSP, since commercial solvers like CPLEX cannot compute the optimal solution for more than fifteen customers due to computation time explosion.In this case, the proposed algorithm includes a sweep strategy as local search method, such that several neighborhoods of the solutions are generated for inserting and swapping orders.A combined drone-vehicle collaborative problem is solved in [29] by, firstly, planning truck and drones' routes independently while accounting for vehicle capabilities and restricted areas, and, secondly, exploiting a genetic search for finding the optimal task scheduling.The research in [30] proposes a genetic algorithm to compute the truck route in between the drone launch sites and a K-means clustering algorithm to compute the optimal launch sites of the UAVs, minimizing the total delivery time and the energy consumed by the vehicles.A genetic algorithm is also proposed in [31] to tackle the formulated mFSTSP that also includes the possibility of carrying Autonomous Transport Vehicles (ATVs) on the truck.Randomly generated instances demonstrate the validity of the approach.A genetic algorithm in combination with a pseudo-random heuristic rule in the crossover phase enhances the exploration capability of the search space in the work proposed in [32].The proposed approach finds both the take-off points and the delivery schedules of the UAVs while minimizing the total completion time using the evolutionary-based strategy.Another genetic algorithm is proposed in [33] to address a scheduling problem of an mFSTSP, minimizing the total distance traveled.In this case, the Clarke and Wright's savings algorithm is used sequentially with the genetic algorithm to enhance the performances.Other genetic-based routing and scheduling solutions can be found in [34][35][36].In [35], a particle swarm algorithm solution is compared to an evolutionary-based solution, showing the superiority of the latter in finding solutions corresponding to a higher level of optimality.The work in [36] implements a multipath genetic algorithm for allocating parcel delivery tasks in emergencies, with some nodes having higher priorities over others.An alternative to evolutionary-based approaches is represented by computationally lighter greedy solutions that, on the other hand, are less likely to produce a more optimal allocation with respect to genetic algorithms.The study in [37] demonstrates the capability of a fast heuristic method to effectively handle an mFSTSP: the methodology combines (i) a multi-armed bandit selection for selecting the truck's customers' partitions from the whole customers' set, (ii) a stochastic local search algorithm that explores the search space, and (iii) a greedy scheduler for estimating the solution quality.The works [38][39][40] show an example of how to address a mFSTSP (in the single drone formulation) and the corresponding heuristic resolution method design in [38], the extension of the formulation to a fleet of multiple drones in [39], and, finally, a further formulation extension for optimizing the drones' speeds within the decision variables domain in [40].These works show how heavily simplifications and assumptions can influence the formulation and the resolution of such problems.A Truck and Drone Routing Algorithm (TDRA) based on Adaptive Large Neighborhood Search (ALNS) is proposed in [41], since the MILP formulation yields to a drastic computation time increase with more than eight customers.It is worth noticing that, beside the optimization methodologies that sub-optimally solve such problems, a way to address the scalability issues of the commercial solvers is to design new mathematical statements of the problem.An example is represented by the study in [25], where a set of mFSTSP formulations are analyzed with the aim of enhancing the size of the largest instances solved in the literature.

Assumptions and Problem Formulation
The single truck-multiple drones routing problem, i.e., the mFSTSP, is the problem of allocating a set of parcel delivery tasks either to a truck or to the multiple UAVs that operate in coordination with it.A parcel delivery task consists of a customer, i.e., a location or node, to be visited exactly once either by the truck or by one of the UAVs of the fleet, which is assumed to be homogenous (all drones have the same characteristics).It is possible that not all of the customers are feasibly serviceable by the UAVs: the set of such customers is defined as the truck only nodes set.The reasons behind the definition of this set may be due to the parcel's weight that exceeds the UAVs' capacity, the difficult landing conditions at the customer's location, the requirement of customer's signature, etc. truck node and drone node denote nodes exclusively served by the truck or a UAV respectively, while drone eligible node is the set of nodes that can be served either by a truck or a UAV.truck route (or tsp tour) denotes, given a set of customers, the minimum cost sequence of nodes visited by the truck.
The delivery process starts with the truck and the drones departing from the depot and ends with the return of all of the vehicles to the depot itself.Each UAV of the fleet can travel toward customers being transported by the truck, conserving its battery life, or can perform a drone sortie, consuming the energy stored in the battery.A drone sortie consists of a set of three nodes from which the UAV is launched, delivers the parcel, and is recovered by the truck.A drone sortie is feasible if the UAV does not exceed its endurance limit; the loaded UAV may be launched from the initial depot or from a customer's node; the recovery node may be another customer's node or the final depot, in the latter case no synchronization between the truck and the UAV is needed since it is assumed that some personnel are always available to recover the drone.Multiple drone sorties can be performed by each UAV in the delivery cycle.drone travel defines a UAV operation moving from one node to another.During the execution of a drone sortie, the truck may go from the launch node to the recovery node or it can accomplish multiple deliveries (accounting for the synchronization constraints that impose on the truck to be present at the recovery node before the battery life of the UAV expires).sub − route denotes the sub-set of the truck route delimited by the node where the drone is launched and the node where the drone is retrieved.The set of all sub − route of a given tsp tour is defined as the truck sub − route.Some truck driver's service times must be considered in order to model a realistic drone-based last-mile delivery service: the service time for, possibly, charging a UAV and loading it with the parcel, a service time to physically deliver the parcel to the customer, and a service time to allow the recovery of a landing UAV.The adopted notation is defined as follows: • Q = {1, 2, . . . ,v}: set of available UAVs for the delivery process.• C = {1, 2, . . . ,c}: set of customers to be served either by the truck or a UAV.• C = {1, 2, . . . ,c } ⊆ C: sub-set of customers that may be served by a UAV.• N = {0, 1, . . . ,c + 1}: set of nodes to be visited exactly once in the delivery process (initial and final node correspond to the same location, i.e., the depot).

•
(i, j, k) q : drone sortie of UAV q, i.e., tuple describing the operation of UAV q being launched at node i, delivering the parcel at node j and being recovered at node k. • P q : set of all possible tuples (i, j, k) q feasibly flown by UAV q.

•
x ij = {0, 1}: decision variable representing whether the truck travels from node i ∈ N 0 to node j ∈ N + (x ij = 1) or not (x ij = 0).• y ijkq = {0, 1}: decision variable representing whether UAV q flies from node i ∈ N 0 to node j ∈ N + and is recovered at node k ∈ {N + : (i, j, k) ∈ P} (y ijkq = 1) or not (y ijkq = 0).• t j : truck arrival time at node j ∈ N + .• t jq : UAV q arrival time at node j ∈ N + .• p ij = {0, 1}: decision variable representing whether customer i ∈ C has been visited before customer j ∈ {C : i = j} by the truck (p ij = 1) or not (x ij = 0).• 1 ≤ u i ≤ c+2: integer position variable of node i ∈ N + in the truck's travel.It follows the mathematical programming-based formulation of the addressed synchronized single truck-multiple UAVs problem, with notation and assumptions stated above.min (t max ) , subject to : (1) 1 Equation ( 1) is the objective function that seeks to minimize the makespan, i.e., the total completion time, of the delivery cycle: it is equivalent to minimize the latest arrival time either of the truck or the latest UAV, i.e., min max t c+1 , t c+1,1 , t c+1,2 , . . ., t c+1,v .This equivalence is justified by Constraints ( 14) and ( 15) that impose a time synchronization among the truck and the UAVs at rendezvous nodes.Constraint (2) imposes the single visit condition for each customer, (3) and (4) impose that the truck starts and ends the delivery cycle at the depot exactly once, and (5) eliminates the truck sub-tours with the bounds of the auxiliary variable u i being set by (28).Constraint (6) imposes that if the truck visits node j, it must depart from node j once the delivery task is performed, ( 7) and ( 8) establish that the UAVs can be launched and retrieved at any node, (9) imposes that if a UAV is launched at node i and retrieved at node k, both nodes i and k must belong to the truck tour, and (10) expresses an analogous condition for the UAVs.Constraint (11) expresses that if a UAV is launched at node i and retrieved at node k, the truck must visit i before k.Constraints ( 12) and ( 13) avoid that either a UAV or the truck are not present at a launch node; analogously, ( 14) and ( 15) avoid the lack of coordination between the truck and the UAVs at a rendezvous node.Constraint (16) regulates the time update of the truck tour if the truck travels from node h ∈ N 0 to node k ∈ N + : the truck arrival time at node k includes the potential launch and recovery time of the UAVs, and the service time to physically deliver the parcel at node k.Constraint (17) regulates the time update for UAV q tour from node i to node j, considering the UAV flight time and the service time to physically deliver the parcel.Constraint (18) includes the travel time of UAV q from node j to node k and the recovery time of the UAV.The endurance limitations of the UAVs are expressed by Constraint (19), which becomes active if the UAVs' travels take place.Constraints ( 20), ( 21), ( 22) establish the correct ordered truck's path.Constraints ( 23) and ( 25) prevent that a UAV is launched from node l before being retrieved at node k.Constraint (24) defines the starting time of the truck and the UAVs.The domain of the decision variables is defined by Constraints from (26) to (31).
The MILP formulated above is inspired by the one formulated in [38].The formulation is modified to include more than one UAV and the service times required to physically deliver the parcels.The problem is NP-hard and a few sub-optimal resolution approaches are presented in the next sections through a set of heuristic algorithms.

Local Search Algorithm
The Local Search Algorithm (LSA) is one of the first and most adopted heuristic algorithms in operation research.It is a heuristic method that evaluates a sub-space of the search space of the problem, an instance of the class of TSPs in this case, until, in general, a sub-optimal solution is found.At each step, the local search algorithm explores the neighboring solutions of a current solution to find a better solution, given an evaluation criterion.It focuses on a single solution and its immediate neighbors, rather than maintaining a population of solutions.It follows the main steps of an LSA:

•
Algorithm Initialization: The algorithm starts the search process from an initial solution.

•
Fitness Evaluation: The fitness or objective function value of the solution is computed; the fitness of the solution represents its quality or suitability.

•
Neighborhood Generation: The algorithm generates solutions neighboring the current solution through small modifications of the current solution.These modifications can include swapping elements, adding and removing elements, or altering values within the encoding of the solution itself.

•
Neighbor Selection and Evaluation: From the set of generated neighboring solutions, one neighbor is selected as the next solution to explore.The selection process can vary based on different criteria, such as selecting the neighbor with the best improvement in fitness value or choosing a neighbor randomly.

•
Solution Update: If the neighbor solution's quality is better than the current solution, the neighbor solution is set as the new current solution.Otherwise, the current solution remains unchanged.

•
Algorithm Termination: The algorithm repeats steps Fitness Evaluation, Neighborhood Generation, Neighbor Selection and Evaluation, and Solution Update until a termination condition is met.The termination condition can be a maximum number of iterations, reaching a satisfactory solution quality, or the lack of solution optimality improvement for a certain number of iterations.
Local search algorithms are particularly useful for solving optimization problems where the objective is to find a solution that satisfies certain constraints or criteria.They are often applied to large search spaces or problems with a high degree of complexity, where exploring the entire solution space is computationally infeasible.However, LSAs may end up being trapped in local optima, which are suboptimal solutions with better optimality level than their immediate neighbors, but with worse optimality level than the global optimum.
Considering the formulated mFSTSP, the proposed LSA starts from the conventional tsp tour, i.e., the routing solution with a single truck and no drones obtained by a TSP solver, then either sequentially removes nodes from the tsp tour and assigns them to one of the UAVs of the fleet or swaps the nodes of the TSP solutions seeking for an improvement of the solution itself.This process is repeated until the LSA does not find any improvement.An example of how the algorithm works is shown in Figure 1.The pseudo-code of the main body of the proposed mFSTSP-related LSA solution is shown in Algorithm 1.The set of drone-eligible nodes C is obtained removing the truck only nodes from the whole node set C. The tsp tour is calculated through a solver that minimizes the total truck travel time given the input C. The computed tsp tour is split into many sub-tours, called sub − route, that includes the truck nodes where the UAVs are launched and recovered.The set of all sub-routes of the truck is the sub − routes vector; it is initialized in Line 3 and progressively modified when UAV assignments are made through the algorithm execution.Analogous considerations hold for i * , j * and k * , which indicate the UAVs' launch, delivery, and recovery nodes.At each iteration and for each UAV, the LSA allocates a delivery task to the UAV with a greater value of max savings.Each UAV is considered independently, as for a single drone FSTSP problem.The savings() function, shown in Algorithm 2, is called every time a drone eligible node j is investigated in order to compute the time reduction related to the removal of node j from the tsp tour.The cost truck() function, shown in Algorithm 3, is called every time a drone eligible node cannot be assigned to a UAV because the UAV is already operating on the considered sub − route.The function evaluates the insertion of node j in a sub − route as a truck node and swaps the truck node if this brings a reduction in the waiting time between the UAV and the truck.The cost uav() function, shown in Algorithm 4, evaluates the potential assignment of a drone eligible node to a uav sortie, if the node j does not belong to a sub − route where UAV travel has already been assigned.The per f orm update() function updates the tsp tour, the truck sub-routes, the UAVs' assignments and the arrival time of the truck at each node of the tsp tour.per f orm update(), shown in Algorithm 5, is called once all drone eligible nodes have been evaluated and if the LSA is able to produce an improvement to the solution.
In order to ensure the synchronization between the UAVs and the truck at the recovery nodes, the time update() function, shown in Algorithm 6, is called within the per f orm update() function.time update() regulates the awaiting of the truck by each UAV and vice versa, updating the truck arrival time vector at each node, denoted by t.For instance, if the UAV's travel time is greater than the corresponding sub − route, the truck should wait for the UAV recovery at a recovery node before serving the next customer.On the other hand, the UAV should hover in place waiting for the truck if the drone travel time is smaller than the corresponding sub − route.The evolutionary-based algorithm is a heuristic algorithm that performs optimization by applying methodologies of exploration of the search space inspired by the evolution mechanisms of nature.This algorithm explores, with a progressive improvement of the optimality level, a sub-set of the solution space of the problem.According to the terminology of this method, a population made by individuals, i.e., either the delivery schedule of the truck or the delivery schedule of the UAVs in this case, gives part of its "genetic" material to the individuals of the next algorithm's iterations in such a way that a cost function's valuable improvement is followed in the next generations.During the reproduction phase, the individuals are modified in order to meet certain fitness criteria: for instance, another heuristic methodology may be adopted to generate a new set of solutions, which partially contain the information, i.e., the genetic material of the solutions of the current iteration and partially derive from the designed further exploration methodology of the search space.It follows the main steps of evolutionary search methods:

•
Algorithm Initialization: The algorithm starts by creating a population of potential solutions to the problem.Each individual of the population is a potential solution to the problem.

•
Fitness Evaluation: Each solution in the population is evaluated and assigned to a fitness value based on how well it solves the problem.The fitness value is represented by the total delivery time of truck and UAVs.

•
Natural Selection: Solutions with higher fitness values are more likely to be selected for the next generation.This process mimics the natural selection of individuals with favorable traits.

•
Individuals Reproduction: The selected solutions are used to create offspring for the next generation.This is done through crossover and mutation operations, which mimic the genetic processes of recombination and unexpected mutation in nature.The crossover phase involves combining genetic material from two parent solutions to create one or more offspring solutions.This is done by selecting a crossover point and exchanging genetic material beyond that point between the parents.The mutation phase introduces small random changes in the genetic material of offspring solutions.It helps explore new regions of the solution space that may lead to better solutions.

•
Individuals Replacement: The offspring solutions replace some of the less fit solutions in the current population.This ensures that the population evolves over time and improves its overall fitness.

•
Algorithm Termination: The algorithm continues to iterate through Fitness Evaluation, Natural Selection, Individuals Reproduction, and Individuals Replacement steps until a termination condition is met.This condition could be a maximum number of generations, a specific fitness threshold, or a predefined computational limit.
By iteratively repeating the steps indicated above, the evolutionary approach explores the solution space, favoring solutions with higher fitness values.Over time, the population tends to converge towards better solutions, but not, in general, to the global optimal solution.The proposed algorithm is a hybrid genetic algorithm that includes optimization and a local search-based heuristic rule during the mutation phase.For this variant of the proposed genetic algorithms, the individuals of the population of each iteration can only correspond to feasible task schedules.The fitness improvement criterion of the solution aims at the minimization of the delivery time of the delivery process.The first proposed variant of the proposed evolutionary-based routing method is the Hybrid Genetic Algorithm-Feasible Variant, shown in Algorithm 7.

Algorithm 7 Main body of the Hybrid Genetic Algorithm-Feasible Variant
Drones 2023, 7, x FOR PEER REVIEW

•
Algorithm Initialization: The algorithm starts by creating a population of solutions to the problem.Each individual of the population is a potential so the problem.

•
Fitness Evaluation: Each solution in the population is evaluated and assig fitness value based on how well it solves the problem.The fitness value is rep by the total delivery time of truck and UAVs.

•
Natural Selection: Solutions with higher fitness values are more likely to be for the next generation.This process mimics the natural selection of individu favorable traits.

•
Individuals Reproduction: The selected solutions are used to create offsprin next generation.This is done through crossover and mutation operation mimic the genetic processes of recombination and unexpected mutation in The crossover phase involves combining genetic material from two parent s to create one or more offspring solutions.This is done by selecting a crossov and exchanging genetic material beyond that point between the pare mutation phase introduces small random changes in the genetic material of o solutions.It helps explore new regions of the solution space that may lead solutions.

•
Individuals Replacement: The offspring solutions replace some of the less fit s in the current population.This ensures that the population evolves over t improves its overall fitness.

•
Algorithm Termination: The algorithm continues to iterate through Fitness Ev Natural Selection, Individuals Reproduction, and Individuals Replacement step termination condition is met.This condition could be a maximum nu generations, a specific fitness threshold, or a predefined computational limi By iteratively repeating the steps indicated above, the evolutionary a explores the solution space, favoring solutions with higher fitness values.Over population tends to converge towards better solutions, but not, in general, to th optimal solution.The proposed algorithm is a hybrid genetic algorithm that optimization and a local search-based heuristic rule during the mutation phase variant of the proposed genetic algorithms, the individuals of the population iteration can only correspond to feasible task schedules.The fitness improvement of the solution aims at the minimization of the delivery time of the delivery pro first proposed variant of the proposed evolutionary-based routing method is th Genetic Algorithm-Feasible Variant, shown in Algorithm 7. Starting from a random population of n individuals, with an individual being defined as a tsp tour without the initial and final depots, the parents selection() function, shown in Algorithm 8, selects the TSP solutions that provide the basis genetic material for the child generation phase.The child generation() function, shown in Algorithm 9, generates new individuals that carry the genetic material of the selected parents with some diversification, to enlarge the solution space explored.The crossover and mutation phases takes place through the split(), local search() and restore() functions that split, modify, and re-assemble individuals in order to generate new genetic material, seeking for an improvement of the solution optimality level with respect to the previous generations.The split() function is reported in Algorithm 10, while its inner functions create uav sorties(), insert in f easible customers(), truck sub − routes update() are reported in Algorithm 11, Algorithm 12, and Algorithm 13 respectively.create uav sorties() assigns all the drones eligible nodes to the UAVs, but if this results in an infeasible assignment, insert in f easible customers() re-inserts the infeasible nodes as truck nodes in the child vector.
The minimum number of customers that can be served by v UAVs, n d , is computed by means of Equation ( 32), with the quantity c−v v+1 being taken from [38] and denoting the minimum number of customers to be served by the truck in presence of v UAVs.
Once all nodes are assigned either to the truck or the UAVs, truck sub − routes update() updates the solution.The population generated by the Hybrid Genetic Algorithm-Feasible Variant considers only feasible solutions, i.e., solutions that do not violate the endurance constraints of the UAVs.When the child has been educated by the local search() function, which is the same as in Algorithm 1, the restore() function, shown in Algorithm 14, generates an individual coherent with the others, since the LSA returns each truck sub − route.The select survivors() function, shown in Algorithm 15, discards the individuals with poor characteristics either in terms of fitness values (the completion time of the delivery cycle) or in terms of diversification from the other individuals.At the end of the whole process, the select f ittest individuals() function selects the best solution according to the minimum makespan criterion.

Algorithm 8 parents selection () function
Drones 2023, 7, x FOR PEER REVIEW Starting from a random population of  individuals, with an individu defined as a   without the initial and final depots, the   function, shown in Algorithm 8, selects the TSP solutions that provide the basis material for the child generation phase.The ℎ () function, sh Algorithm 9, generates new individuals that carry the genetic material of the parents with some diversification, to enlarge the solution space explored.The c and mutation phases takes place through the () ,  ℎ() and  functions that split, modify, and re-assemble individuals in order to generate new material, seeking for an improvement of the solution optimality level with respe previous generations.The () function is reported in Algorithm 10, while functions   () ,   () ,   () are reported in Algorithm 11, Algorithm 12, and Algor respectively.  () assigns all the    to the but if this results in an infeasible assignment,   () r the infeasible nodes as   in the ℎ vector.
The minimum number of customers that can be served by  UAVs,  , is co by means of Equation ( 32), with the quantity being taken from [38]

Hybrid Genetic Algorithm-Infeasible Variant
The second proposed variant of the evolutionary-based solution approac formulated mFSTSP is the Hybrid Genetic Algorithm-Infeasible Variant, sh Algorithm 16.The main difference with the Hybrid Genetic Algorithm-Feasible V that the requirement of feasible solution output of the () function, sh Algorithm 17, does not hold, and the definition of the additional () f shown in Algorithm 18, is needed.The second variant considers infeasible ind with the aim of diversifying the exploration of the search space and avoiding in local minima.

Hybrid Genetic Algorithm-Infeasible Variant
The second proposed variant of the evolutionary-based solution approach of the formulated mFSTSP is the Hybrid Genetic Algorithm -Infeasible Variant, shown in Algorithm 16.The main difference with the Hybrid Genetic Algorithm-Feasible Variant is that the requirement of feasible solution output of the split() function, shown in Algorithm 17, does not hold, and the definition of the additional repair() function, shown in Algorithm 18, is needed.The second variant considers infeasible individuals with the aim of diversifying the exploration of the search space and avoiding incurring local minima.increment variable reduces the number of nodes randomly assigned to the UAVs.This contains the number of infeasible solutions in the customers node set if the node locations are too far to be served by the UAVs.The create uav sortie() function assigns all the customers to the uav sortie independently on their feasibility.

Ad-hoc Method
The alternative solution method Ad-hoc Method, shown in Algorithm 19, is a greedy algorithm that differs from the other optimization algorithms that iteratively try to improve a set of sub-optimal solutions, and generates a unique solution with a systematic procedure.A greedy algorithm is a simple, intuitive, and often efficient algorithmic approach that makes locally optimal choices at each step, aiming at finding an optimal solution.At each step of the algorithm, a well-defined strategy builds the best possible solution at the current step, without considering the overall consequences of that choice.It follows a general outline of how a greedy algorithm works:

•
Algorithm Initialization: The algorithm starts with an empty or partial solution.

•
Greedy Strategy: At each step, the algorithm makes a locally optimal choice based on a certain rule.It selects the most favorable available at the current stage of the process, without considering its impact on future steps.

•
Solution Feasibility and Update: The algorithm checks if the selected choice is feasible or satisfies certain constraints.If the choice violates any constraints, it may be discarded, and an alternative option may be considered.The selected choice is added to the current solution or modifies the partial solution.

•
Algorithm Termination: The algorithm iterates through Greedy Strategy and Solution Feasibility and Update steps until a termination condition is met.This condition could be reaching a final solution, satisfying a specific objective, or exhausting the available choices.
Greedy algorithms do not guarantee an optimal solution for all problems, as in this case.Therefore, the correctness and optimality of a greedy algorithm depend on the specific problem being solved, the methodology exploited by the greedy algorithm itself and the properties of the problem instance.Greedy algorithms are often used when finding an optimal solution for a problem requires a lot of computational resources or is not feasible within a reasonable time frame.Additionally, greedy algorithms are commonly used as building blocks or components within more complex algorithmic approaches.
The greedy approach proposed in this work for the formulated mFSTSP is inspired by the work in [38], which also formulated a greedy strategy for a mFSTSP.The proposed greedy solution consists of a poor computationally expensive ad-hoc resolution method of the formulated mFSTSP problem that computes a feasible schedule solution according to a specific methodology of synchronization between the UAVs' delivery schedules and the truck's delivery schedule.served by the truck if  UAVs are available.Maximizing the number of custom served by the UAVs aims at reducing the total completion time as much as poss at the same time,  may be increased for feasibility reasons if the maximum nu customers cannot be assigned to the UAVs.' and   are initialized at ev start and the   () function, shown in Algorithm 20, out shortest possible   given  customers.Firstly, the  −   inserted in the   , since they cannot be assigned to any UAV; secon   are inserted in the   from the  .
The  () function decouples the nodes to be served by t and the nodes to be served by the UAVs.Then, the proposed algorithm tentativel LTL (Lower Truck Limit) denotes the minimum number of customers that can be served by the truck if v UAVs are available.Maximizing the number of customers to be served by the UAVs aims at reducing the total completion time as much as possible and, at the same time, LTL may be increased for feasibility reasons if the maximum number of customers cannot be assigned to the UAVs.C and tsp tour are initialized at every loop start and the cost truck calculation() function, shown in Algorithm 20, outputs the shortest possible tsp tour given LTL customers.Firstly, the truck − only nodes are inserted in the tsp tour, since they cannot be assigned to any UAV; secondly, the customers nodes are inserted in the tsp tour from the truck input.
The customer division() function decouples the nodes to be served by the truck and the nodes to be served by the UAVs.Then, the proposed algorithm tentatively creates a sub-optimal mFSTSP solution according to the available tsp tour and uav customers, and repeats the process increasing LTL if the solution is not feasible or if a reduction in the makespan can be made, even if the solution is feasible.

Simulation Settings
The proposed greedy algorithms are implemented by means of the Py programming language.In particular, the DEAP library is exploited for management of the genetic algorithm variables.The settings of the evol algorithms' parameters are as follows:  = 1000,  = 20,  = 1,  = 15, The parameter settings for the fleet of UAVs are reported in Table 1.Note tha customers nodes are considered eligible for the UAVs.As far as the optimizatio process' makespan, this choice allows us to highlight the benefits of drones' cooperation with a truck.

Simulation Settings
The proposed greedy algorithms are implemented by means of the Python 3.9 programming language.In particular, the DEAP library is exploited for ease of management of the genetic algorithm variables.The settings of the evolutionary algorithms' parameters are as follows: γ = 1000, n = 20, ω = 1, η = 15, µ = 25.The parameter settings for the fleet of UAVs are reported in Table 1.Note that all the customers nodes are considered eligible for the UAVs.As far as the optimization of the process' makespan, this choice allows us to highlight the benefits of drones' usage in cooperation with a truck.The delivery operational area of either the truck or the drones is defined as a rectangular portion of the city of Turin, Italy.The fleet is assumed to be composed of v = 2 UAVs.δ t is set to 200 seconds.The truck-only TSP solutions are obtained with the Google OR tools solver [42], with the parameter local search metaheuristic set as guided local search and the parameter time limit set to 30 s.The capability of the proposed approaches of managing the formulated mFSTSP is tested considering two scenarios and the results are corroborated by Monte Carlo simulations.Scenario A consists of c = 10 customers to be served (plus the depot, denoted with 0, from which the truck departs and returns to end the delivery cycle), randomly sampled in the operational area.Table 2 shows the customers' coordinates.Such coordinates are given to the Bing Maps API REST Services [43] in order to obtain the truck travel times and the drones travel times, encoded as cost matrixes.Scenario B is defined analogously, this time considering c = 20 customers to be served in a larger rectangular area defined in the city of Rome, Italy.This is done to test the proposed solutions with a more challenging scenario, with a higher number of customers and greater travel distances.Table 3 shows the randomly sampled customers' coordinates, the corresponding cost matrixes are derived with the same procedure adopted for Scenario A.

Results and Discussion
Considering either Scenario A or Scenario B, the improvement that is obtained in terms of total completion time adopting the UAVs in combination with the truck with respect to the conventional truck-only TSP solution is shown in Table 4. Considering a single iteration of the proposed algorithms, Table 4 also shows the level of optimality in terms of the makespan of the delivery cycle that can be achieved by each routing method for the two scenarios of references.The results show that the enhanced solution space exploration capability of the genetic algorithms results in the highest savings over the TSP solutions among the proposed algorithms, regardless of the instance size.This is also confirmed in the works [31][32][33][34] where evolutionary-based methods outperform most of the alternative approaches.
Nevertheless, the computational complexity of the evolutionary approach is not justified for small size instances.This is confirmed by the results of Scenario A: the greedy method generates a delivery schedule with a slightly lower optimality level (~3%) with respect to the evolutionary methods, but with a noticeable gain of computational time.This is no longer the case when a bigger instance size is considered, as shown by the results of Scenario B. The differences between the two variants of the evolutionary approach emerge in Scenario B, as the feasible variant is more likely to find a better optimal solution when the maximum number of iterations is limited.When the search space is broader, the genetic material of feasible solutions that is passed from one population to the next one represents a sub-optimal search direction of the solution space.If the size of the problem instance increases, this is most likely to produce a better optimal solution with respect to the infeasible variant.For the sake of clarity, a visual representation of the solutions obtained with the TSP solver and the proposed routing algorithms is given in Figures 2-6 for Scenario A, starting from the conventional TSP solution up to the Ad-hoc Method solution.Each of the mentioned figures shows the truck and drones' travels on the map of the portion of the city of Turin, Italy, and the corresponding Gantt Chart with the timed schedules.Considering either Scenario A or Scenario B, one hundred Monte Carlo simulations are performed for each scenario in order to corroborate the capability of the proposed evolutionary algorithms to handle the formulated mFSTSP.Table 5 shows the results.The infeasible variant of the evolutionary-based task scheduling method has enhanced exploration capability of the solution space and it is less likely to stop at local minima too early.This is demonstrated with Scenario A, where both the minimum and the maximum makespan results are smaller than the ones related to the feasible variant.Nevertheless, with a more complex scenario and a bigger operational area, i.e., as in Scenario B, the optimality level of the solution produced by the Hybrid Genetic Algorithm-Infeasible Variant is worse.This is due to the fact that, when the search space is much broader as for this scenario, the constraint of considering only feasible schedules during the crossover phase restricts the investigated solutions to the neighborhood of the feasible solutions and this is more likely to generate a more optimal schedule, given a "reasonable" maximum number of iterations of the genetic algorithms.As expected, and accordingly to the search strategies of the two evolutionary-based methods, the variability of the results is smaller for the feasible variant than for the infeasible variant.By paying the price of heavier computational loads, both the variants of the evolutionary-based methods ensure a greater solution's optimality level than either the greedy Local Search Algorithm or the Adhoc Method, for both Scenario A and Scenario B.  Considering either Scenario A or Scenario B, one hundred Monte Carlo simulations are performed for each scenario in order to corroborate the capability of the proposed evolutionary algorithms to handle the formulated mFSTSP.Table 5 shows the results.The infeasible variant of the evolutionary-based task scheduling method has enhanced exploration capability of the solution space and it is less likely to stop at local minima too early.This is demonstrated with Scenario A, where both the minimum and the maximum makespan results are smaller than the ones related to the feasible variant.Nevertheless, with a more complex scenario and a bigger operational area, i.e., as in Scenario B, the optimality level of the solution produced by the Hybrid Genetic Algorithm-Infeasible Variant is worse.This is due to the fact that, when the search space is much broader as for this scenario, the constraint of considering only feasible schedules during the crossover phase restricts the investigated solutions to the neighborhood of the feasible solutions and this is more likely to generate a more optimal schedule, given a "reasonable" maximum number of iterations of the genetic algorithms.As expected, and accordingly to the search strategies of the two evolutionary-based methods, the variability of the results is smaller for the feasible variant than for the infeasible variant.By paying the price of heavier computational loads, both the variants of the evolutionary-based methods ensure a greater solution's optimality level than either the greedy Local Search Algorithm or the Ad-hoc Method, for both Scenario A and Scenario B. As far as the problem addressed in this work is concerned, the choice of the most suitable method between the feasible and the infeasible variants depends on the instance of the problem itself and on the available computational resources.Considering both the inherent features of the proposed routing algorithms and the simulation results, the summarizing considerations of this work are reported in Table 6.

Conclusions
In this paper, a few heuristic methods were designed, implemented, and compared for a last-mile delivery (LMD) Travelling Salesman Problem with a truck and multiple drones.Firstly, a mathematical formulation of the problem was proposed considering the case of a truck that cooperates synchronously with a fleet of Unmanned Aerial Vehicles.The problem formulation accounts for most of the key representative elements of an LMD service, such as the delivery due dates of the parcels, the service time required to process the orders, the battery capacities of the UAVs, etc.The goal was to find a combined truck-UAVs feasible schedule to minimize the total completion time of the delivery process while respecting the constraints related to either the delivery tasks or the vehicles' capacities.Secondly, given the addressed NP-hard problem, we designed a local search solution approach, two variants of an evolutionary-based solution and an ad-hoc heuristic solution specifically designed for addressing the formulated problem.The latter was conceptually different from the other proposed methods as it generated a solution systematically, considering the features of the specific problem.The evolutionary-based strategy was split into two variants depending on the feasibility feature of the solutions that were permutated at each execution step.The proposed solutions were evaluated over two realistic urban air mobility scenarios.The capability of the methods to produce an optimal schedule was corroborated through Monte Carlo simulation campaigns.
Future research directions will include the extension of the problem formulation to multiple trucks and heterogenous drones and the design of a dynamic task allocation strategy in order to account for environmental uncertainties and order modifications during task execution.The design and implementation of other routing methods may also be investigated and compared to the algorithms proposed in this work.

Algorithm 1 1 .
Main body of the Local Search Algorithm Algorithm Main body of the Local Search Algorithm.Initialize:  , ,   , ,  ,  ,  , ;   = 0;  = \  ;   −  = [ [ ] as many  are available ];  * = [ [] as many  are available ] ;  * = [ [] as many  are available ] ;  * = [ [] as many  are available ];     = [ [] as many  are available ]; while no improvement is produced do for  ∊  do for  ∊  do (,, ); for  −  ∊   − [] do if  −  is already assigned to a UAV do  (,,,   − ,  − ); else  (,,, − ); end end end if   > 0 do  ( * [],  * [], * [],  [],   − ,  ); either sequentially removes nodes from the   and assigns them to one of the UAVs of the fleet or swaps the nodes of the TSP solutions seeking for an improvement of the solution itself.This process is repeated until the LSA does not find any improvement.An example of how the algorithm works is shown in Figure 1.The pseudo-code of the main body of the proposed mFSTSP-related LSA solution is shown in Algorithm 1.

Figure 1 .
Figure 1.Representation of the functioning logic of the LSA with a truck and one drone.Grey nodes represent the depot, red nodes represent truck-only nodes and green nodes represent drone-eligible nodes.(a) The TSP tour computed by the TSP solver; (b) 1st iteration of the LSA: node 4 is removed from the TSP tour and assigned to the UAV, which is launched at node 2 and retrieved at node 3; (c) 2nd iteration of the LSA: node 5 is removed from the TSP tour and assigned to the UAV, which is launched at node 3 and retrieved at the depot.

Drones 2023, 7 ,Figure 2 .Figure 3 .
Figure 2. TSP solution representations of Scenario A. (a) TSP tour reported on the map of the operational area, with the paths of the truck in blue; (b) Gantt chart of the TSP solution.The numbers from 1 to 10 in the map and the Gant chart represent the Customer Node IDs.

Figure 2 .Figure 2 .Figure 3 .
Figure 2. TSP solution representations of Scenario A. (a) TSP tour reported on the map of the operational area, with the paths of the truck in blue; (b) Gantt chart of the TSP solution.The numbers from 1 to 10 in the map and the Gant chart represent the Customer Node IDs.

Figure 3 .
Figure 3. Local Search Algorithm solution representations of Scenario A. (a) mFSTSP tour reported on the map of the operational area, with the paths of the truck (blue), drone 1 (orange) and drone 2 (red); (b) Gantt chart o f the mFSTSP solution.The numbers from 1 to 10 in the map and the Gant chart represent the Customer Node IDs.

Figure 4 .Figure 5 .
Figure 4. Hybrid Genetic Algorithm-Feasible Variant solution representations of Scenario A. (a) mFSTSP tour reported on the map of the operational area, with the paths of the truck (blue), drone 1 (orange) and drone 2 (red); (b) Gantt chart of the mFSTSP solution.The numbers from 1 to 10 in the map and the Gant chart represent the Customer Node IDs.

Figure 4 .
Figure 4. Hybrid Genetic Algorithm-Feasible Variant solution representations of Scenario A. (a) mFSTSP tour reported on the map of the operational area, with the paths of the truck (blue), drone 1 (orange) and drone 2 (red); (b) Gantt chart of the mFSTSP solution.The numbers from 1 to 10 in the map and the Gant chart represent the Customer Node IDs.

Figure 4 .Figure 5 .
Figure 4. Hybrid Genetic Algorithm-Feasible Variant solution representations of Scenario A. (a) mFSTSP tour reported on the map of the operational area, with the paths of the truck (blue), drone 1 (orange) and drone 2 (red); (b) Gantt chart of the mFSTSP solution.The numbers from 1 to 10 in the map and the Gant chart represent the Customer Node IDs.

Figure 5 .Figure 6 .
Figure 5. Hybrid Genetic Algorithm-Infeasible Variant solution representations of Scenario A. (a) mF-STSP tour reported on the map of the operational area, with the paths of the truck (blue), drone 1 (orange) and drone 2 (red); (b) Gantt chart of the mFSTSP solution.The numbers from 1 to 10 in the map and the Gant chart represent the Customer Node IDs.

Figure 6 .
Figure 6.Ad-hoc Method solution representations of Scenario A. (a) mFSTSP tour reported on the map of the operational area, with the paths of the truck (blue), drone 1 (orange) and drone 2 (red); (b) Gantt chart of the mFSTSP solution.The numbers from 1 to 10 in the map and the Gant chart represent the Customer Node IDs.

Table 2 .
[44]omer nodes' coordinates for Scenario A (10 customers and 1 depot), generated randomly via[44]in a rectangular map comprehending the Turin urban area.

Table 3 .
[44]omer nodes' coordinates for Scenario B (20 customers and 1 depot), generated randomly via[44]in a rectangular map comprehending the Rome urban area.

Table 4 .
Makespan and savings over truck-only TSP of each proposed solution for Scenario A and Scenario B.

Table 5 .
Monte Carlo simulation results: total completion time of the delivery process.

Table 5 .
Monte Carlo simulation results: total completion time of the delivery process.

Table 6 .
Summarizing consideration about the proposed routing algorithms.