Optimization of a Mixed Fleet of Aerial Drones for Medical Supplies: A Case Study of Blood Delivery Logistics

Abstract

Aerial drones have emerged as an innovative solution for faster transportation of time-sensitive items (e.g., emergency medical supplies), potentially reducing the transmission of contagious diseases and enhancing healthcare availability through contactless autonomous delivery. We study fleet sizing and efficient scheduling of a mixed fleet of drones for delivering time-sensitive medical items having distinct release and due times to minimize the required fleet size and fleet composition, the required number of additional batteries, and the total energy consumption. We continuously track the remaining battery energy of drones to determine the optimal timing for battery replacement, rather than replacing the battery at each node. Using actual drone flight test data, we employed a machine learning (ML) method to estimate the energy consumption of different drone types during flight segments for different operating parameters. We present a novel mixed-integer programming model to efficiently formulate the problem that integrates the estimated energy consumption functions from ML. We propose a new greedy heuristic (GH) algorithm and a customized genetic algorithm (GA) for solving large-scale instances of this problem faster. Results demonstrate that the GH algorithm is substantially faster than the accelerated CPLEX and the GA, while sacrificing the solution quality by a small amount. Results based on an actual blood sample delivery case study from Pendleton, Oregon, United States, show that using a mixed fleet of drones reduces the total cost and total energy consumption up to 18.18% and 28.7%, respectively, compared to using a homogeneous fleet.

Motivation

The benefits of aerial delivery include but are not limited to overcoming the limitations in ground transportation, minimizing the possibility of cross-infection through contactless delivery, and reducing delivery times. Recent studies have increasingly acknowledged the potential of aerial drones in facilitating medical assistance and enhancing emergency response durations ( 1 – 3 ). Studies have demonstrated that drones can improve humanitarian aid supply chain amidst the COVID-19 pandemic, specifically for transporting emergency medical items, such as blood samples, vaccines, drugs, pharmaceuticals, laboratory specimens, medical equipment, and test kits ( 4 , 5 ). Examples of successfully using drone technology in real-world medical deliveries include rural Virginia ( 6 ) and Ghana ( 7 ). In addition, the United States Postal Service and Zipline have collaborated to assess the feasibility of delivering medications, blood, and vaccines in Rwanda ( 8 ). Other developing nations have undertaken similar initiatives (e.g., Mahimbo [ 9 ] and Mis [ 10 ]).

However, drone-based delivery is constrained by its battery capacity, necessitating efficient scheduling and routing decisions, despite its significant potential. Current studies on drone deployment and operation typically assume either a fixed flight duration ( 11 ), energy consumption based solely on drone speed and package weight ( 12 ), and/or assume drones can make several deliveries in a single flight. These simplifying assumptions overlook drones’ energy consumption, potentially leading to impractical routing/scheduling, in-flight failures causing severe accidents, and increased customer dissatisfaction. However, in practice, not only the drone speed and package weight, but also flight path, travel distance, and components of flight duration (e.g., ascend, descend, and hover) affect the energy consumption of drones ( 13 ), often ignored in current studies.

Furthermore, existing studies often assume drone batteries are charged/replaced at each node in the delivery network without considering the drone’s remaining energy and the ability to perform additional deliveries, leading to an unnecessarily large number of additional batteries, and/or drones. Moreover, recharging drone batteries during delivery operations is generally impractical for time-sensitive medical logistics having strict release and due times due to the significantly long recharging time. For instance, charging the battery of a DJI drone takes 51–115 min ( 14 ), whereas it takes approximately 5 min to replace the battery of a DJI Matrice 600 Pro ( 13 ). Due to much shorter battery replacement time, drones’ batteries can be replaced with fully charged batteries during the logistics operations and recharged overnight at a sustainable rate when the drones are not being used. In practice, drones are not operated at night for delivery logistics. Incorporating battery recharging during delivery operations of time-sensitive medical items can lead to severe consequences (e.g., patient deaths) due to delivery delay in critical health-related logistics (e.g., blood products and emergency medical supplies) and/or a considerably larger number of drones to satisfy tight delivery time windows due to reduced drone utilization caused by long recharging time. Moreover, frequently recharging the battery of the drones during the logistics operations degrades the battery life and makes the batteries reach their end of life faster ( 15 ). In contrast, an endogenous battery replacement strategy, where the drone battery is replaced only when it is needed accounting for the delivery sequence, continuously keeping track of the remaining battery energy, and the minimum required battery energy, reduces unnecessary battery replacements and the required number of additional batteries. This strategy enhances drones’ utilization and ensures cost-efficient, timely deliveries in satisfying the strict time windows (defined by release and due times) of time-sensitive medical items.

Another common simplification is assuming drones fly in straight paths, despite Federal Aviation Administration (FAA)-designated “No Drone Zones” (e.g., schools, airports, large buildings, and restricted private property) requiring detours via waypoints, increasing energy consumption and delivery time. Additionally, most existing studies consider a homogeneous fleet of drones with identical characteristics. However, in medical delivery logistics, demands are usually generated from geographically distributed healthcare business locations (i.e., delivery/pickup locations for drones) that are at different distances from the depot, and have different package weights (e.g., blood sample versus defibrillator), release (i.e., order ready) times, and due times. A mixed fleet comprising different drone types with distinct characteristics—cost, battery capacity, power consumption, speed, and package weight carrying capacity—is more efficient to address this variation in delivery requirements, as different drone types are suitable for different delivery tasks, such as different delivery distances, package weights, and release times and time windows, due to their unique characteristics. For instance, rotary drones are generally more cost- and energy-efficient than fixed-wing vertical-takeoff-and-landing (VTOL) drones for delivering lighter packages to shorter distances, whereas VTOL drones are larger with larger battery capacity and offer a more cost-efficient solution than rotary drones for long-distance deliveries with heavier packages ( 13 ). Using a mixed fleet allows leveraging the benefits of different drone types to make the logistics operation more cost- and energy-efficient than a homogeneous fleet. Therefore, the varying characteristics and requirements (i.e., package weight, distance, and release and due times) of geographically distributed healthcare business locations (i.e., delivery/pickup locations for drones) necessitate leveraging a mixed fleet of drones with distinct characteristics to make the time-sensitive medical items delivery logistics cost- and energy-efficient.

Therefore, it is critical to make up and route a mixed fleet of drones and schedule their battery replacement efficiently accounting for practical factors to minimize fleet size and determine the best mixed fleet composition, and to minimize the required number of additional batteries and total energy consumption while delivering the time-sensitive medical items maintaining the specified release and due times. However, routing and scheduling a mixed fleet of drones with endogenous battery replacement increases the computational complexity of the decision-making process due to a larger number of alternative routing and scheduling solutions. This necessitates developing efficient mathematical model and solution algorithms to quickly and accurately solve this complex logistics problem, especially for large-scale problem instances.

Related Literature

Scientific research on using drones in different areas of logistics is experiencing rapid growth, as presented in recent studies (e.g., Moshref-Javadi and Winkenbach [ 16 ] and Macrina et al. [ 17 ]). Some researchers have centered on providing an overview of operations research techniques utilized in addressing last-mile delivery challenges through the use of drones ( 18 ), while others have directed their attention toward reviewing the optimization of drone applications in civilian contexts, such as the delivery of medical supplies and disaster management ( 19 ).

The current optimization modeling domain extensively explores optimizing medical supply transportation using drones. Rabta et al. ( 20 ) presented a theoretical framework for humanitarian logistics in a network consisting of multiple nodes and a single depot. The proposed model aimed to minimize drone flight distance while considering payload and energy limits. Wen et al. ( 21 ) presented a model for the capacitated vehicle routing problem (VRP) that utilizes drones in delivering emergency blood supplies. The model aims to minimize the costs associated with distance and flight time. Chowdhury et al. ( 22 , 23 ) studied drone routing optimization considering battery recharging strategy. The authors constructed mathematical models to determine the most advantageous placement of distribution centers and the expenses incurred in post-disaster evaluations linked to these variables. Gentili et al. ( 24 ), Kim et al. ( 25 ), and Shavarani et al. ( 26 ) studied the facility location problem with drones. Shavarani et al. ( 26 ) proposed an optimization model with two objectives, namely the minimization of the overall costs of facility location and costs related to drone operation and the reduction of the number of customers left uncovered. Kim et al. ( 25 ) examined the stochastic nature of flight distance variability resulting from fluctuations in battery consumption within a modeling framework. The optimization problem in Gentili et al. ( 24 ) focused on an optimization problem that aimed to minimize the overall disutility value of emergency medical supplies accounting for their perishability. Their study assumes that each platform has a drone serving a customer on each flight. The optimization of drone battery endurance in medical supply logistics has been the subject of research by Dhote and Limbourg ( 27 ), Ghelichi et al. ( 28 ), and Macias et al. ( 29 ). These studies have concentrated on the joint optimization of a charging station or tactical center location selection and drone route, with consideration given to the battery life of the drones. However, most of these studies have yet to look at the realistic energy consumption of drones using actual flight test data and the effect of different flight paths on energy consumption. Moreover, these studies did not explicitly model battery replacement strategies. Furthermore, they overlooked the effect of using a mixed fleet of drones to explore the efficiency of drone-based delivery systems in healthcare logistics and emergency medical supply.

Existing literature on drone-based delivery has identified various logistics models, including pure-play drone-based models, unsynchronized multi-modal models, synchronized multi-modal models, and re-supply multi-modal models ( 16 ). Among them, pure-play drone-based models are relevant to the issue of drone-only delivery, which involves the direct transportation of packages from a warehouse to customers utilizing drones. The alternative models specifically address the issue of truck-drone coordination in the context of delivery services (e.g., Jeong et al. [ 30 ], Kuo et al. [ 31 ], and Murray and Chu [ 32 ]). Most studies focusing on drone-only delivery problems generally operate under the assumption that multiple drones are available within the warehouse and that each drone can serve one or more nodes in a single trip. Dorling et al. ( 12 ) presented a pair of models for the VRP in the dynamic delivery planning context. The first model aims to minimize overall operating expenses by satisfying delivery time limitations, while the second model aims to optimize delivery time under budgetary constraints. The model employed by the authors utilized a linear approximation function for computing power consumption. Notably, this function did not consider the distinct impact of various flight components of drones on energy consumption. However, Cheng et al. ( 33 ) extended this research by formulating power consumption as a nonlinear function of payload and travel time while also considering time windows. Although the logical and subgradient cut methods were employed to handle the nonlinear power function, the authors did not verify the computed energy consumption through actual drone flight test data. Choi and Schonfeld ( 34 ) conducted a study into an autonomous drone delivery system in which all customers’ requirements were identical. The proposed model aims to optimize drone fleet size considering the interplay between battery capacity, package weight, and flight range. However, the authors did not compute the energy consumption of drones accounting for practical factors, such as flight path, and different flight components.

Most of the previous studies suffer from the limiting assumption that a drone battery can complete each trip, delivering multiple packages, without considering practical battery replacement. Additionally, previous studies mostly assume that the drone fleet is primarily homogeneous, ignoring the potential energy-efficient delivery by employing a mixed fleet. In a recent study, Wang et al. ( 35 ) considered using a mixed fleet of drones in a drone-based delivery system. However, the authors did not provide insights into the system’s energy efficiency concerning using a mixed fleet of drones, nor did they consider the fleet sizing and fleet composition. Furthermore, many research studies have not accounted for the realistic battery replacement, where the drone battery replacement decision is made endogenously based on the drone routing and scheduling decisions. This limiting assumption may lead to an increased number of unnecessary battery replacements, the need for additional batteries, and charging infrastructure as well as less drone utilization and a larger required fleet size. Recently, only Bhuiyan et al. ( 13 ) have incorporated the drone battery replacement process as a means to minimize the frequency of battery replacements. But the authors considered battery replacement at the depot only and did not consider battery replacement at the delivery locations along with the depot. Additionally, their drone scheduling model is developed for a homogeneous fleet. Furthermore, the authors compared the effects of no-fly zones on energy consumption, fleet size, and battery replacement by considering flying over the road network, which is not always necessary. The authors did not consider the realistic flight path where drones mostly follow a straight path while taking detours following some waypoints to avoid no-fly zones. Completely flying over road networks to avoid no-fly zones increases the energy consumption, fleet size, and required battery replacements substantially compared to taking necessary detours via some waypoints following straight paths to avoid these no-fly zones. Only Jeong et al. ( 30 ) integrated this realistic flight path considering circular no-fly zones into the drone’s flight path. However, they did not account for the realistic battery replacements, nor did they model the problem considering a mixed fleet of drones.

Our study is relevant to the drone routing/scheduling for pickup and delivery services, which is beneficial for delivering medical supplies, especially during a pandemic to avoid human involvement. As our study incorporates a machine learning (ML)-based drone energy consumption in the mixed-fleet drone scheduling model for delivering/picking up packages with release and due dates, it can easily address the drone scheduling for simultaneous pickup and delivery. Existing research studies addressed the drone routing/scheduling for pickup and delivery with drone-only (e.g., Wikarek et al. [ 36 ] and Gómez-Lagos et al. [ 37 ]) and drone-truck fleet (e.g., Meng et al. [ 38 ]). Ham ( 39 ) applied a constraint programming approach to expand the parallel drone scheduling traveling salesman problem by incorporating drop and pickup synchronization, where after completing a delivery, a drone can either return to the depot to deliver the next customer order or directly fly to the next customer for pickup. Pachayappan and Sudhakar ( 40 ) studied a drone pickup and delivery problem using a docking station for communication and resource support. Wikarek et al. ( 36 ) studied another variant of the pickup and delivery problem with drones where a drone can both deliver and pick up a package from the same customer. Gómez-Lagos et al. ( 37 ) studied another variant of the pickup and delivery drone routing problem where each drone leaves a depot to pick up a package from a facility and then deliver to a customer location. All of these studies ignored the drone reuse in a multi-trip setting, leading to solutions requiring an unnecessarily large number of drones in delivering time-sensitive packages with release and due dates. Recently, Meng et al. ( 38 ) studied a multi-trip drone routing problem with trucks for pickup and delivery services. A summary of the most relevant existing studies on drone-based healthcare logistics is presented in Table 1.

Table 1.

A Summary of the Existing Studies on Drone Deliveries for Healthcare Logistics

Ref.	Mixed fleet	Drone energy consumption				NFZ	ADD	EBR
		Package weight	Distance	Drone type	Different flight segments
( 1 )	×	×	×	×	×	×	×	×
( 2 )	×	✓	✓	×	×	×	✓	×
( 12 )	×	✓	✓	×	✓	×	✓	×
( 20 )	×	✓	✓	×	✓	×	×	×
( 21 )	×	✓	✓	×	×	×	×	×
( 22 )	×	✓	✓	×	✓	×	×	×
( 23 )	✓	✓	✓	✓	✓	×	×	×
( 24 )	×	×	✓	×	×	×	×	×
( 25 )	×	×	✓	×	×	×	×	×
( 27 )	×	✓	✓	×	×	✓	×	×
( 28 )	×	✓	✓	×	×	×	×	×
( 41 )	×	✓	✓	×	×	×	×	×
( 42 )	×	×	×	×	×	×	×	×
( 43 )	×	×	×	×	×	×	×	×
( 44 )	×	×	×	×	×	×	×	×
This paper	✓	✓	✓	✓	✓	✓	✓	✓

Note: Ref. = reference; NFZ = no-fly zone; ADD = actual drone data; EBR = endogenous battery replacement.

“✓” and “×” $\mathrm{denote}$ whether the corresponding paper considered the associated aspect or not, respectively.

In summary, despite the valuable insights from the existing literature, several significant research gaps persist. First, the most important research gap in the existing literature is the lack of routing and scheduling a mixed fleet of aerial drones with endogenous battery replacement for delivering time-sensitive healthcare items having distinct release times and time windows. This is particularly critical because a homogeneous fleet of drones is not cost- and energy-efficient to deliver time-sensitive items across geographically distributed locations with different characteristics (i.e., delivery distance, package weight, and release time and time window). Moreover, considering endogenous battery replacement of drones, where instead of replacing drone battery at each node in the delivery network, battery replacement is decided endogenously in the model accounting for the routing/scheduling decisions, continuously keeping track of the remaining battery energy, and the minimum required battery energy, is overlooked in the existing studies. Second, simultaneous consideration of practical factors in computing drone energy consumption, including package weight, drone type, distance, and different flight segments (e.g., hover, ascend, descend, and forward flight) based on actual drone flight test data as well as accounting for no-fly zones for the flight path of drones in the transportation network are overlooked in the existing studies. Third, most recent studies lack an efficient mixed-integer linear programming (MILP) model and fast algorithms for efficiently computing routing and scheduling solutions for a mixed drone fleet incorporating the above practical aspects with real data on drones.

Contributions

In this research, we aim to (1) develop a mathematical optimization model and algorithms for business owners to efficiently make up and operate a mixed fleet of drones for the direct delivery of time-sensitive medical items; and (2) provide insights into how various business and drone operating conditions affect the required fleet composition (i.e., the required number of drones of each type in the fleet), as well as the cost and energy consumption for healthcare logistics and emergency medical supplies. Specifically, the contributions of this paper are as follows.

First, we study a novel problem of leveraging a mixed fleet of aerial drones with distinct characteristics—cost, battery capacity, speed, package weight carrying capacity, and power consumption—and endogenous battery replacement for collecting/delivering time-sensitive medical items having distinct release times and due times from/to geographically distributed healthcare business locations (e.g., clinics). A mixed fleet of drones offers a more cost- and energy-efficient logistics system compared to a homogeneous drone fleet, particularly in real-life logistics operations, where demand points are at different distances from the depot, have different package weights, and require the packages to be picked up and delivered at different times. Moreover, we endogenously decide the drone battery replacement in the problem accounting for the routing/scheduling decisions, continuously keeping track of the remaining battery energy, and the minimum required battery energy.

Second, we simultaneously consider different practical factors, including package weight, drone type, distance, and different flight segments in computing the energy consumption of drones for different flight paths—straight path, no-fly zone, and road network.

Third, we develop a novel MILP model strengthened by a valid inequality that incorporates these practical aspects to determine the optimal routing and scheduling of the mixed fleet that minimizes the required fleet size, the required number of additional batteries, and the total energy consumption of drones. By explicitly considering minimizing the cost of required additional batteries in the objective function, our model tries to minimize the required number of additional batteries in the logistics operation, which is critical for budget-constrained applications, especially when specialized drone batteries have considerably high costs. We also propose a new greedy heuristic (GH) algorithm and a problem-specific customized genetic algorithm (GA) to solve the problem for large-scale instances faster. Moreover, we utilized an ML method to estimate drone energy consumption for different package weights in different flight segments using actual drone flight data. These estimated values of the energy consumption of drones are then used in the MILP model and solution algorithms as parameters.

Fourth, we present numerical results and new managerial insights based on an actual blood sample delivery case study from the Interpath Laboratory, Inc.—a healthcare logistics company—located in Pendleton, Oregon, United States, into the effect of various practical drone operating parameters, such as package weight, flight path, the minimum required battery energy, fleet type, and time window (defined by the release and due times) on the required fleet composition (i.e., required number of drones of each type), required fleet size, required number of additional batteries, and total energy consumption.

Problem Description

This study examines the challenge faced by a healthcare business entity, such as a medical research laboratory, in optimally using a mixed (i.e., heterogeneous) fleet of drones to collect/deliver time-sensitive medical items (e.g., blood samples, blood products, test kit) with distinct release and due times from/to other business organizations, such as hospitals or clinics. Figure 1 shows a visualization of the mixed drone fleet routing and scheduling to serve geographically distributed healthcare business locations. The healthcare business owner seeks to determine the best composition of the mixed drone fleet (i.e., number of drones of each type), fleet size, and additional batteries in conducting its logistics operations for time-sensitive medical items. This business-to-business drone delivery problem studied in this paper is applicable to diverse industries, such as the transportation of medical equipment, emergency relief goods delivery, and small parts delivery in automobile manufacturing. We assume the business owner owns the mixed fleet of drones and the medical laboratory serves as the depot for the drone fleet. We refer to this research facility as the central depot throughout this paper. The drone fleet consists of different types of drones having different physical characteristics, such as speed, package weight carrying capacity, battery capacity, and power consumption during the flight segments. Other healthcare business locations (e.g., hospital, clinic) center around this depot. Each business location has battery swapping and charging infrastructure. Hospitals or clinics place orders to the depot at different times throughout a planning horizon. We assume that the service provider knows these order placement times, locations, and package weights before the planning horizon starts. The service provider seeks to schedule the drone fleet based on the order placement time, location, and the package weight of each order. We consider that each order is associated with a pickup time window, defined by the order release time (i.e., package ready time) and the due time (i.e., timestamp before which the package should be picked up or delivered). Different orders can have different release and due times and thus different time windows; for instance, some medical deliveries may have a narrow time window of only 30 min, whereas others may allow a more flexible pickup time window of several hours.

Figure 1.

Visualization of the routing and scheduling a mixed fleet of drones to collect/deliver time-sensitive medical items from/to geographically distributed healthcare business locations.

In this healthcare delivery logistics problem, each drone directly flies to the pickup/delivery location and then returns to the depot before flying to the next location. We assume that the drones have a fully charged battery at the beginning of the planning horizon. However, the battery’s energy gradually decreases as drones perform pickup and delivery operations. The energy consumption rate depends on various operating factors, such as drone speed, package weight, travel distance, duration of flight segments, and flight path. To ensure the safe operation of the drones, before leaving for the next destination (next node in the logistics network), we must ensure that the drone battery has at least a minimum required amount of energy remaining after arriving at that location. Therefore, before leaving the current location (either the depot or other business location) for the next destination, we need to compute the potential remaining energy in the drone battery after arriving at the following location. If the potential remaining energy is less than the minimum required energy, we must replace the current battery with a fully charged one before departing from the current location. This flexibility of replacing batteries in each node in the logistics network significantly extends the delivery range of a drone type, which is particularly critical in the healthcare logistics. However, starting with a fully charged battery from a node (either from the depot or a pickup/delivery location), if a drone’s remaining energy falls below the minimum required energy after reaching the next destination node, then this drone type cannot serve the destination node starting from the origin node. In other words, this origin-destination arc is not feasible in the logistics network for this drone type that depends on the drone speed, package weight, flight path, battery capacity, and minimum required battery energy. We refer to these destination nodes as outside the delivery range for the drone type. Additionally, each drone type has a maximum package weight-carrying capacity, and a drone type cannot serve a location for which the package weight exceeds this maximum package weight-carrying capacity. Accounting for all these factors, the service provider’s goal is to optimally schedule the mixed fleet of drones to deliver/pick up the customer orders maintaining their release and due times with the best composition of the mixed fleet, as well as find the minimum required fleet size, additional batteries, and total cost and energy consumption.

To improve clarity and readability, we list main assumptions of our problem as follows.

• A mixed fleet of drones with distinct characteristics, including cost, battery capacity, power consumption, package weight-carrying capacity, and speed, is stationed at a central depot (e.g., medical laboratory) with a fixed location.

• A set of healthcare business organizations (e.g., hospitals) with fixed locations as delivery/pickup locations that are geographically located around the central depot, where drones should collect/deliver time-sensitive medical items.

• Each healthcare business location places an order to the depot, where each order is associated with a distinct package weight, package release time, and pickup/delivery due time.

• Each drone carries a single package on a trip to a delivery location. Moreover, each drone performs multiple direct delivery trips between the central depot and different healthcare business locations, where a trip starts with a drone picking up a package from the central depot, delivering it to a healthcare business location, and then returning to the depot to conclude this trip.

• Each drone starts its route from the depot with a fully charged battery and replaces its battery before departing for the next location in an as-needed fashion accounting for the routing/scheduling decisions, continuously keeping track of the remaining battery energy, and the minimum required battery energy.

• All nodes in the network (i.e., the central depot and healthcare business locations) have available fully charged batteries to replace a depleted battery as needed to support safe and efficient operation of drones. A battery replacement for each drone type takes a fixed amount of time.

• To ensure practicality, some drone types cannot carry/deliver certain packages that exceed the drone’s maximum package weight carrying capacity, and/or the energy consumption of the drone exceeds its battery capacity in delivering the packages to the corresponding healthcare business locations.

Mathematical Formulations

To solve the mixed fleet of drones scheduling problem with release and due times and endogenous battery replacement discussed in the section “Problem description” we formulate the problem as an MILP model. In this section, we present the estimation of drones’ energy consumption using an ML method, the MILP model, and the valid inequality to improve the computational efficiency of the MILP model. We list the necessary sets, parameters, and variables that support the mathematical formulations as follows.

Set Description

$\mathcal{I}$  Set of delivery/pickup locations $i$

${\mathcal{I}}^{\prime }$  Set of delivery/pickup locations $i$ including dummy source $0$ and sink $S$

$\mathcal{D}$  Set of drone types indexed by $d$

${\mathcal{I}}_d$  Set of delivery/pickup locations $i\in {\mathcal{I}}^{\prime }$ that drone type $d\in \mathcal{D}$ can serve

Parameter Description

${\pi }_d$ Amortized hourly cost of a drone of type $d$

${\pi }_{\mathrm{bat}}^d$ Amortized hourly cost of a battery of drone type $d$

${\pi }_{\ell }$ Hourly cost of each drone operator

$v_d$ Speed of drone type $d$

${\pi }_E$ Cost per unit of energy

$n_d$ Number of drones a drone operator can operate simultaneously

${\pi }_d^M$ Maintenance cost of drone type $d$

$D_i$ Distance of delivery/pickup location $i$ from the depot

$t_{\ell }^d$ Loading time for a drone of type $d$

$t_u^d$ Unloading time for a drone of type $d$

$t_{0i}^d$ Total time required for a drone of type $d$ to arrive at location $i$ from the depot

$t_{i0}^d$ Total time required for a drone of type $d$ to return from location $i$ to the depot

$t_{\mathrm{bat}}^d$ Time required to replace the battery of drone type $d$

$t_{\mathrm{as}}^d$ Time required to ascend at flying height from the ground for a drone of type $d$

$t_{\mathrm{fl}}^d$ Forward flight time for a drone of type $d$ to fly from depot to location $i$

$t_h^d$ Required hover duration before descending for a drone of type $d$

$t_{\mathrm{ds}}^d$ Time required to descend to the ground from flying height for a drone of type $d$

$t_{\mathrm{rn}}$ Time required to hover while taking a turn to change direction

$n_{\mathrm{rn}}$ Required number of turns to arrive location $i$ from depot with no-fly zone detouring

$b_d^{\max }$ Initial energy in the battery of drone type $d$

$b_d^{\min }$ Minimum remaining energy required in the battery of drone type $d$

$E_{0i}^d$ Energy consumption by a drone of type $d$ to arrive at node $i$ from the depot

$E_{i0}^d$ Energy consumption by a drone of type $d$ to return from node $i$ to the depot

$E_{\mathrm{ij}}^d$ Energy consumption by a drone to fly from node $i$ to $j$ , that is, $E_{\mathrm{ij}}^d=E_{i0}^d+E_{0j}^d$

$w_d^{\max }$ Maximum package weight carrying capacity of a drone of type $d$

$w_i$ Package weight of node $i$

Variable Description

$x_{\mathrm{ij}}^d$ 1, if locations $i$ is served immediately before $j$ by drone type $d$ , 0 otherwise

$y_{0i}^d$ 1, if battery of drone type $d$ is replaced after arriving at node $i$ , $0$ otherwise

$y_{i0}^d$ 1, if battery of drone type $d$ is replaced after returning to the depot from node $i$ , $0$ otherwise

$q_{d0i}$ Remaining energy in the battery of drone type $d$ after arriving at node $i$

$q_{\mathrm{di}0}$ Remaining energy in the battery of drone type $d$ after returning to the depot from node $i$

$q_{d0i}^{\prime }$ Potential remaining battery energy of drone type $d$ after arriving at node $i$

$q_{\mathrm{di}0}^{\prime }$ Potential remaining battery energy of drone type $d$ after returning to the depot from node $i$

$p_i^d$ Timing of when a drone of type $d$ picks-up/delivers the package for node $i$

Estimation of Energy Consumption of Drones

As realistic and accurate estimation of energy consumption is crucial for an efficient and robust drone-based delivery system, we used a supervised ML method to estimate drones’ energy consumption based on our actual drone flight data. The section “Case study” provides a detailed overview of the related data used to calculate the energy consumption of drones considered in this study. In this section, we discuss the methods used to estimate the total energy consumption in an entire flight of a drone type. Each flight of a drone consists of the following flight segments—ascend, descend, forward flight, and hover. Drones consume different power and energy in different flight segments. The power and energy consumption of a drone vary with package weight, speed, and other operational parameters. We performed different drone flight tests for each type of drone with different package weights, speeds, turn angle and radius, duration and patterns of flight segments (i.e., ascend, descend, and hover), and other operating parameters, among which the individual effect of package weight on power consumption is the largest. As there are numerous different package weights in the delivery data, we need to estimate the energy consumption as a function of package weight for the ease of use in the drone scheduling MILP model. We trained our ML model based on the power consumption data from our drone flight tests to estimate the power consumption in different flight segments for different package weights. The power consumption data for each package weight used in the ML model account for the specific values of other parameters, such as speed, duration and patterns of flight segments (i.e., ascend, descend, and hover), turn angle and radius. We found that the linear regression model is the most appropriate choice to estimate power consumption in each flight segment for different package weights. The adjusted $R$ -squared values for the model are greater than 98%, confirming a linear relationship exists between power consumption and package weight. Figure 2 shows the actual and estimated power consumption in each flight segment of a hexacopter (i.e., DJI Matrice 600 Pro) drone type. We denote the fitted functions estimating the power consumption for a given package weight $m$ during ascend, descend, forward flight, and hover, as $P_{\mathrm{as}}^d\left(m\right)$ , $P_{\mathrm{ds}}^d\left(m\right)$ , $P_{\mathrm{fl}}^d\left(m\right)$ , $P_h^d\left(m\right)$ , respectively. Therefore, the total time and energy consumption of a drone type $d$ in a single trip from the depot to a particular location $i$ while flying over a straight path can be calculated as follows.

Figure 2.

Estimation of power consumption: (a) ascend, (b) descend, (c) forward flight, (d) hover.

Total time required to arrive at $i$ from the depot, $t_{0i}^d=t_{\mathrm{as}}^d+t_{\mathrm{fl}}^d+t_h^d+t_{\mathrm{ds}}^d$ . As the drone speed is the same in the forward and return trips, the total time required to return to the depot from location $i$ , $t_{i0}^d=t_{0i}^d$ . Total energy consumption to arrive empty (i.e., $m=0$ ) at location $i$ from the depot, $E_{0i}^d=P_{\mathrm{as}}^d\left(0\right)\times t_{\mathrm{as}}^d+P_{\mathrm{fl}}^d\left(0\right)\times t_{\mathrm{fl}}^d\times D_i+P_h^d\left(0\right)\times t_h^d+P_{\mathrm{ds}}^d\left(0\right)\times t_{\mathrm{ds}}^d$ . Similarly, the total energy consumption to return to the depot from location $i$ with a package of weight $m$ , $E_{i0}^d=P_{\mathrm{as}}^d\left(m\right)\times t_{\mathrm{as}}^d+P_{\mathrm{fl}}^d\left(m\right)\times t_{\mathrm{fl}}^d*D_i+P_h^d\left(m\right)\times t_h^d+P_{\mathrm{ds}}^d\left(m\right)\times t_{\mathrm{ds}}^d$ .

In the presence of no-fly zones in the flight path, a drone needs to detour following some waypoints, resulting the drone to make turns in changing its flight directions. While making these turns to change directions, drones usually have to hover for a short duration. These turns and hover due to avoiding no-fly zones result in an increased drone travel time and energy consumption compared to flying in a straight path between a source-destination pair. We include this additional travel time and energy consumption due to the no-fly zone detouring with the travel time and energy consumption of flying in a straight path to ensure an accurate energy consumption modeling. Therefore, when a drone needs to avoid no-fly zones following some waypoints to serve a particular location $i$ , the total time and energy consumption are calculated as follows:

Total time required to arrive at $i$ from the depot avoiding no-fly zones, $t_{0i}^d=t_{\mathrm{as}}^d+t_{\mathrm{fl}}^d+t_h^d+t_{\mathrm{ds}}^d+t_{\mathrm{rn}}\times n_{\mathrm{rn}}$ . Total time required to return to the depot from location $i$ , $t_{i0}^d=t_{0i}$ . Total energy consumption to arrive empty (i.e., $m=0$ ) at location $i$ from the depot avoiding no-fly zones, $E_{0i}^d=P_{\mathrm{as}}^d\left(0\right)\times t_{\mathrm{as}}^d+P_{\mathrm{fl}}^d\left(0\right)\times t_{\mathrm{fl}}^d\times D_i+P_h^d\left(0\right)\times t_h^d+P_{\mathrm{ds}}^d\left(0\right)\times t_{\mathrm{ds}}^d+P_h^d\left(0\right)\times t_{\mathrm{rn}}^d\times n_{\mathrm{rn}}^d$ . Similarly, the total energy consumption to return to the depot from location $i$ with a package of weight $m$ , $E_{i0}^d=P_{\mathrm{as}}^d\left(m\right)\times t_{\mathrm{as}}^d+P_{\mathrm{fl}}^d\left(m\right)\times t_{\mathrm{fl}}^d\times D_i+P_h^d\left(m\right)\times t_h^d+P_{\mathrm{ds}}^d\left(m\right)\times t_{\mathrm{ds}}^d+P_h^d\left(m\right)\times t_{\mathrm{rn}}^d\times n_{\mathrm{rn}}^d$ .

MILP Model

As mentioned in the section Problem description, in this mixed-fleet drone scheduling problem, each drone directly flies from the depot to a pickup/delivery location and returns to the depot before flying to the next location, meaning each drone makes multiple trips from the depot over the planning horizon. To efficiently solve this multi-trip mixed-fleet drone scheduling problem with release and due times (MTMFDSP-RD) by an exact method, we transformed and modeled this problem as a VRP. To accomplish this, we introduced a dummy source ( $0$ ) and a sink ( $S$ ), representing the depot, to the set of delivery locations $\mathcal{I}$ , resulting in a new set ${\mathcal{I}}^{\prime }$ . Accounting for the attributes of the pickup/delivery tasks and the corresponding energy consumption by the available drone types, we compute the set of pickup/delivery locations, ${\mathcal{I}}_d$ , drone type $d\in \mathcal{D}$ can serve. To model this multi-trip mixed-fleet drone scheduling problem as VRP, we define a “route” for a drone that represents the sequence of locations visited by the drone over the planning horizon. For instance, consider in the entire planning horizon, a drone of type $d$ visits locations $1$ , $3$ , and $4$ sequentially such that the original drone paths are as follows: drone first visits location $1$ from the depot, then returns to depot before visiting location $3$ , and after returning from location $3$ to the depot, the drone visits location $4$ and returns to the depot. In our new “route” concept, the route $0-1-3-4-S$ (shown in the example network in Figure 3 as VRP transformation path) represents the sequence of locations visited by the drone. This route $0-1-3-4-S$ also represents a drone.

Figure 3.

Example network representation of vehicle routing problem (VRP) transformation.

Leveraging this “route” concept and the sets ${\mathcal{I}}_d$ , we can model the MTMFDSP-RD by defining variables for drone type $d$ , avoiding the need to define a variable for each individual drone. We define a variable $x_{\mathrm{ij}}^d$ to represent whether a drone of type $d$ visits location $j$ immediately after returning to the depot from the location $i$ or not. Note that although $d$ represents the drone type, each individual drone follows a distinct route. The corresponding parameter to this variable, $E_{\mathrm{ij}}^d$ , denotes the energy consumed by the drone type $d$ to return to the depot from the location $i$ and then arrive at the location $j$ to perform the next delivery/pickup operation; that is, $E_{\mathrm{ij}}^d=E_{i0}^d+E_{0j}^d$ . The parameter $E_{\mathrm{ij}}^d$ for each drone type is estimated using the ML method (as described in the section “Estimation of energy consumption of drones”) to use in the MILP model. The variable $x_{0j}^d$ represents an outgoing arc from the dummy source (0) and the beginning of a route (or the sequence of deliveries) by a drone of type $d$ . The number of outgoing arcs from the dummy source for each drone type $d$ represents the total number of drones of the drone type $d$ used. Therefore, adding a large penalty to the corresponding energy consumption parameter, $E_{0j}^d$ , we can minimize the number of outgoing arcs (i.e., routes) from the source and thus minimize the required fleet size, which is the number of drones of each type used in the healthcare delivery logistics.

We present the MILP model (1) for the MTMFDSP-RD below.

$$\begin{array}{c}\multicolumn{1}{c}{\min \{\sum _{d\in \mathcal{D}}\sum _{i\in {\mathcal{I}}^{\prime }}\sum _{j\in {\mathcal{J}}^{\prime }}{\pi }_E{E}_{\mathrm{ij}}^d{x}_{\mathrm{ij}}^d+\sum _{d\in \mathcal{D}}\sum _{j\in {\mathcal{J}}^{\prime }}({\pi }_d+\frac{{\pi }_{\ell }}{n_d}+{\pi }_d^M)x_{0j}^d}\\ \multicolumn{1}{c}{+\sum _{d\in \mathcal{D}}\sum _{i\in {\mathcal{J}}^{\prime }}{\pi }_{\mathrm{bat}}^d(y_{0i}^d+y_{i0}^d)\}}\end{array}$$

(1a)

$$s.t.\sum _{d\in \mathcal{D}}\sum _{i\in {\mathcal{I}}^{\prime }}{x}_{\mathrm{ij}}^d=1,\forall j\in \mathcal{J}$$

(1b)

$$\sum _{d\in \mathcal{D}}\sum _{j\in {\mathcal{J}}^{\prime }}{x}_{\mathrm{ij}}^d=1,\forall i\in \mathcal{J}$$

(1c)

$$\sum _{d\in \mathcal{D}}\sum _{i\in {\mathcal{J}}^{\prime }}{x}_{\mathrm{ij}}^d=\sum _{d\in \mathcal{D}}\sum _{i\in {\mathcal{J}}^{\prime }}{x}_{\mathrm{ji}}^d,\forall j\in \mathcal{J}$$

(1d)

$$\sum _{j\in \mathcal{J}}{x}_{0j}^d=\sum _{j\in \mathcal{J}}{x}_{\mathrm{jS}}^d,\forall d\in \mathcal{D}$$

(1e)

$$\sum _{i\in {\mathcal{I}}^{\prime }}{E}_{0j}^d{x}_{\mathrm{ij}}^d\le (b_d^{\max }-b_d^{\min }),\forall d\in \mathcal{D},j\in {\mathcal{J}}_d$$

(1f)

$$w_j\sum _{i\in \mathcal{J}}{x}_{\mathrm{ij}}^d\le w_d^{\max },\forall d\in \mathcal{D},j\in {\mathcal{J}}_d$$

(1g)

$$\begin{array}{c}\multicolumn{1}{c}{p_j^d\ge p_i^d+(t_{i0}^d+t_{\ell }^d+t_u^d+t_{0j}^d)+t_{\mathrm{bat}}^d(y_{0i}^d+y_{i0}^d)}\\ \multicolumn{1}{c}{-M(1-x_{\mathrm{ij}}^d),\forall d\in \mathcal{D};i,j\in {\mathcal{J}}^{\prime },i\ne j}\end{array}$$

(1h)

$$e_i\le p_i^d\le {\ell }_i,\forall d\in \mathcal{D},i\in \mathcal{J}$$

(1i)

$$q_{d0j}^{\prime }\le b_d^0-E_{0j}^d+M(1-x_{0j}^d),\forall d\in \mathcal{D},j\in {\mathcal{J}}_d\setminus \{0\}$$

(1j)

$$q_{\mathrm{dj}0}^{\prime }\le q_{d0j}-E_{j0}^d,\forall d\in \mathcal{D},j\in {\mathcal{J}}_d\setminus \{0\}$$

(1k)

$$q_{d0j}^{\prime }\le q_{\mathrm{di}0}-E_{0j}^d+M(1-x_{\mathrm{ij}}^d),\forall d\in \mathcal{D},i,j\in {\mathcal{J}}_d,i\ne j$$

(1l)

$$q_{\mathrm{dj}0}^{\prime }-b_d^{\min }\ge M_{\ell }^d{y}_{0j}^d,\forall d\in \mathcal{D},j\in {\mathcal{J}}_d\setminus \{0\}$$

(1m)

$$q_{\mathrm{dj}0}^{\prime }-b_d^{\min }\le M_u^d(1-y_{0j}^d),\forall d\in \mathcal{D},j\in {\mathcal{J}}_d\setminus \{0\}$$

(1n)

$$\begin{array}{c}\multicolumn{1}{c}{q_{d0j}^{\prime }-b_d^{\min }\ge M_{\ell }^d{y}_{i0}^d+M(x_{\mathrm{ij}}^d-1),\forall d\in \mathcal{D},}\\ \multicolumn{1}{c}{i\in {\mathcal{J}}^{\prime },j\in {\mathcal{J}}_d\setminus \{0\},i\ne j}\end{array}$$

(1o)

$$\begin{array}{c}\multicolumn{1}{c}{q_{d0j}^{\prime }-b_d^{\min }\le M_u^d(1-y_{i0}^d)+M(1-x_{\mathrm{ij}}^d),}\\ \multicolumn{1}{c}{\forall d\in \mathcal{D},i\in {\mathcal{J}}^{\prime },j\in {\mathcal{J}}_d\setminus \{0\},i\ne j}\end{array}$$

(1p)

$$q_{\mathrm{dj}0}\le b_d^{\max }-E_{j0}^d+M_u(1-y_{0j}^d),\forall d\in \mathcal{D},j\in {\mathcal{J}}_d\setminus \{0\}$$

(1q)

$$q_{\mathrm{dj}0}\le q_{d0j}-E_{j0}^d+M_u^d{y}_{0j}^d,\forall d\in \mathcal{D},j\in {\mathcal{J}}_d\setminus \{0\}$$

(1r)

$$\begin{array}{c}\multicolumn{1}{c}{q_{d0j}\le b_d^{\max }-E_{0j}^d+M_u^d(1-y_{i0}^d)+M(1-x_{\mathrm{ij}}^d),}\\ \multicolumn{1}{c}{\forall d\in \mathcal{D},i\in {\mathcal{I}}_d,j\in {\mathcal{J}}_d\setminus \{0\},i\ne j}\end{array}$$

(1s)

$$\begin{array}{c}\multicolumn{1}{c}{q_{d0j}\le q_{\mathrm{di}0}-E_{0j}^d+M_u^d{y}_{i0}^d+M(1-x_{\mathrm{ij}}^d),}\\ \multicolumn{1}{c}{\forall d\in \mathcal{D},i\in {\mathcal{J}}_d\setminus \{0\},j\in {\mathcal{J}}_d\setminus \{0\},i\ne j}\end{array}$$

(1t)

$$\begin{array}{c}\multicolumn{1}{c}{q_{d0j}\le b_d^{\max }-E_{0j}^d+M_u^d{y}_{00}^d+M(1-x_{0j}^d),}\\ \multicolumn{1}{c}{\forall d\in \mathcal{D},j\in {\mathcal{J}}_d\setminus \{0\}}\end{array}$$

(1u)

$$x_{\mathrm{ij}}^d\in \{0,1\},\forall d\in \mathcal{D},i,j\in {\mathcal{J}}^{\prime },i\ne j$$

(1v)

$$y_{0i}^d,y_{i0}^d\in \{0,1\},\forall d\in \mathcal{D},i\in {\mathcal{J}}^{\prime }$$

(1w)

$$q_{0j}^d,q_{j0}^d\ge 0,\forall d\in \mathcal{D},j\in {\mathcal{J}}_d$$

(1x)

$$q_{d0j}^{\prime },q_{\mathrm{dj}0}^{\prime }\ge -b_d^{\max },\forall d\in \mathcal{D},j\in {\mathcal{J}}_d$$

(1y)

$$p_i^d\ge 0,\forall d\in \mathcal{D},i\in {\mathcal{J}}^{\prime }$$

(1z)

The objective function (1a) seeks to minimize the total cost due to energy consumption, drone fleet, and additional batteries in delivering/picking up all the orders satisfying the release and due times. As the cost of drone fleet and additional batteries are functions of the required number of drones and battery replacements, the optimal solution of this model determines the required fleet composition, required fleet size, and the required number of additional batteries to deliver/pick up all orders satisfying their release and due times.

Constraints (1b) and (1c) ensure that each customer order is delivered/picked up once by only one drone of any type. Constraints (1d) are the flow balance constraints for the delivery/pickup locations, ensuring that if a drone visits a location $i$ , it must return to the depot from that location to visit the next location in the delivery sequence until the drone completes its route (i.e., the delivery sequence). Constraints (1e) are the flow balance constraints for the dummy source–sink pair.

Delivery to a location outside the range of a drone is not possible, as stated in the section Problem description. Constraints (1f) enforce the delivery range restriction for each drone type. In addition, constraints (1g) ensure the package weight carrying capacity restriction of each drone type. Constraints (1h) ensure the accurate computation of the pickup/delivery time of each location in the delivery sequence (i.e., route) of a drone. Constraints (1i) enforce the time window for each order, where $e_i$ and ${\ell }_i$ are the release and due times, respectively, for location $i$ .

To ensure safe operation of the drones, there should be a minimum required energy remaining in the drone battery during the drone flight. As in this problem, drone battery can be replaced both at the depot and at the delivery/pickup locations, we need to compute how much energy would remain in the drone’s battery after arriving at the next destination, which can be either the depot or the delivery/pickup location. If the potential remaining energy after arriving to the next destination is below the minimum required energy, the drone battery must be replaced before flying for this next destination. Constraints (1j) compute the potential remaining energy in the drones’ batteries after arriving at the first delivery/pickup location from the depot in each delivery sequence (i.e., route). Each constraint (1l) computes the potential remaining energy in the drone’s battery after arriving at each delivery/pickup location except for the first location in each delivery sequence. Constraints (1k) compute the potential remaining energy after returning to the depot from delivery/pickup location $j$ . If the potential remaining energy after returning from the delivery/pickup location $j$ is less than the minimum required energy (i.e., $q_{\mathrm{dj}0}^{\prime }<b_d^{\min }$ ), we need to replace the current battery with a fully charged one at location $j$ before the drone leaves this location to return to the depot (i.e., $y_{0j}^d=1$ ). Constraints (1m) represent this condition. But if the potential remaining energy after returning from the delivery/pickup location $j$ is higher than the minimum required energy (i.e., $q_{\mathrm{dj}0}^{\prime }\ge b_d^{\min }$ ), then the battery must not be replaced before returning to the depot, and thus, the drone should continue its operation with the remaining energy $\left(q_{\mathrm{dj}0}\right)$ to return to the depot. Constraints (1n) enforce this condition. Similarly, if the potential remaining energy after arriving to delivery/pickup location $j$ from the depot drops below the minimum required energy (i.e., $i.e.,q_{d0j}^{\prime }<b_d^{\min }$ ), the drone battery must be replaced before departing the depot (constraints (1o)). But, the battery must not be replaced before departing the depot for location $j$ if $q_{d0j}^{\prime }\ge b_d^{\min }$ , which is represented by constraints (1p).

Depending on whether the battery is replaced or not before returning to the depot from location $j$ , we need to compute the actual energy remaining in the drone’s battery after returning to the depot from location $j$ . Constraints (1q) ensure that if the drone’s battery is replaced at location $j$ before returning to the depot (i.e., $y_{0j}^d=1$ ), the actual energy remaining after returning to the depot is computed by subtracting the energy consumption due to return-trip to the depot from the initial battery energy, that is, $q_{\mathrm{dj}0}\le b_d^{\max }-E_{j0}^d$ . But, if the drone’s battery is not replaced before leaving location $j$ (i.e., $y_{0j}^d=0$ ), the actual energy remaining after returning to depot from location $j$ is computed based on the remaining energy ( $q_{d0j}$ ) after arriving previously at this location $j$ from the depot as $q_{\mathrm{dj}0}\le q_{d0j}-E_{j0}^d$ , which is represented by constraints (1r). Similarly, constraints (1s) ensure that if the drone’s battery is replaced at the depot after arriving from the previous location $i$ , that is, $y_{i0}^d=1$ (in other words, before leaving the depot for the next location $j$ ), the actual remaining energy after arriving at location $j$ is computed by subtracting the energy consumption in arriving to $j$ (i.e., $E_{0j}^d$ ) from the initial battery energy as $q_{d0j}\le b_d^{\max }-E_{0j}^d.$ In contrast, if the drone’s battery is not replaced before leaving the depot for the next location $j$ (i.e., $y_{i0}^d=0$ ), the actual energy remaining after arriving at location $j$ is computed based on the remaining energy after arriving to the depot from the previous location $i$ (i.e., $q_{\mathrm{di}0}$ ) as $q_{d0j}\le q_{\mathrm{di}0}-E_{0j}^d$ , which is represented by constraints (1t). Each constraint (1u) computes the actual energy remaining after arriving at the first location in each delivery sequence or route. Constraints (1v)–(1z) represent the binary nature and sign restrictions of the decision variables.

Valid Inequality

In the proposed MILP model (1), the number of $x_{\mathrm{ij}}^d$ variables increases by $\left|{\mathcal{I}}^0\right|\times \left|\mathcal{D}\right|$ for each additional delivery/pickup location and adds many constraints to the model, increasing computational complexity substantially. To address this challenge and reduce the computational burden, we introduce a problem-specific valid inequality that eliminates unnecessary $x_{\mathrm{ij}}^d$ variables from the solution space, leading to an improvement in computational efficiency. A drone cannot pick up/deliver the package for location $j$ immediately after serving location $i$ , if the due time for $j$ ( ${\ell }_j$ ) is earlier than the earliest time the drone can arrive at location $j$ after serving location $i$ . Each constraint (2) removes the $x_{\mathrm{ij}}^d$ variables from the candidate solutions for which this condition holds.

$$\begin{array}{c}\multicolumn{1}{c}{\mathrm{if}{e}_i+t_{\ell }^d+t_u^d+t_{i0}^d+t_{0j}^d>{\ell }_j}\\ \multicolumn{1}{c}{\mathrm{then}{x}_{\mathrm{ij}}^d=0\forall i,j\in \mathcal{J}}\end{array}$$

(2)

Solution Algorithms

Due to the non deterministic polynomial time (NP)-hardness of the classic capacitated VRP with time windows, drone-based VRP with additional constraints for battery range is also NP-hard ( 45 ). Our proposed MILP model (1) for a mixed fleet of drones can be easily reduced to capacitated VRP with time windows by relaxing the battery range limitation and assuming a homogeneous fleet, and thus proving it to be NP-hard. Due to having an exponential runtime, our proposed MILP model (1) is computationally challenging to solve for large-scale delivery networks. Therefore, we propose a new GH algorithm and a customized GA to solve the MTMFDSP-RD with larger instances much faster than the off-the-shelf solvers, such as CPLEX.

Greedy Heuristic Algorithm

This subsection presents a GH algorithm for solving the proposed MTMFDSP-RD problem for larger instances. In this algorithm, we used the proposed “route” concept discussed earlier in the section “MILP model,” where a route represents the sequence of delivery/pickup locations served by a drone. Additionally, a route is equivalent to a drone of a particular type. The complete pseudocode of the proposed heuristic algorithm is presented in Algorithm 1. The algorithm starts by sorting the delivery/pickup locations in an ascending order of their package ready time (i.e., the earliest possible pickup time or the release time). Accounting for the energy consumption and package weight carrying capacity of the drone types, we compute a set ${\mathcal{I}}_d$ for each drone type $d\in \mathcal{D}$ that contains the delivery/pickup locations drone type $d$ can serve.

Using these sets ${\mathcal{I}}_d$ , the algorithm finds the drone types that are eligible (based on the drone’s range and package weight-carrying capacity) to serve each delivery/pickup location $i\in \mathcal{I}$ . For each location $i\in \mathcal{I}$ , the algorithm iteratively evaluates the drone types eligible to serve $i$ based on the additional cost due to serving $i$ and selects the best drone type that results in the minimum cost consisting of the cost of energy, battery, and drones. In evaluating each eligible drone type $d$ for serving $i$ , the algorithm iteratively checks whether location $i$ can be served by any of the existing drones of type $d$ (i.e., whether location $i$ can be added to any of the existing routes of that drone type $d$ ) accounting for the release time and required time window of customer order for location $i$ , with the possibility of replacing the drone’s battery. The algorithm evaluates the cost of all the existing routes of drone type $d$ and selects the route with the minimum cost. If multiple routes have the same minimum cost, we select the route with the minimum cost and minimum remaining energy after serving location $i$ . However, if none of the existing routes can serve the location $i$ , we assign a new drone of type $d$ (i.e., create a new route starting with location $i$ from the depot) to serve location $i$ . After finding the best routes (i.e., drones) for each eligible drone type, the algorithm selects the best route of the drone type $d^{*}$ having the minimum cost to serve the location $i$ and then updates the route and its remaining energy accordingly. The algorithm terminates once it completes evaluation for all the locations and returns the required fleet size and composition, total energy consumption, and total number of battery replacements. The procedures of making drone battery replacement decisions and checking the required time window satisfaction in each route are presented as separate functions in Algorithm 2 that are used in Algorithm 1.

Algorithm 1.

Greedy Heuristic Algorithm

Input Input Data Output ${\mathcal{L}}_{\mathrm{rd}},{\mathrm{nBat}}_d^{*}$ , ObjectiveValue 1: procedureGreedyHeuristic(Input Data) 2: Initialize, list of routes, ${\mathcal{L}}_{\mathrm{rd}}=[[0,],]$ , list of drone remaining energies, ${\mathcal{L}}_{\mathrm{Ed}}=[{\mathrm{Ch}}_d^0,]$ , and number of battery replacements $\mathrm{nBa}{t}_d=0$ for each drone type $d\in \mathcal{D}$ 3: Sorted $\mathcal{I}\leftarrow$ sort delivery/pickup locations in ascending order of package ready time 4: ${\mathcal{I}}_d\leftarrow$ set of delivery/pickup locations that can be served by drone type $d$ 5: ObjectiveValue $\leftarrow 0$ 6: for each delivery, $i$ , in the sorted $\mathcal{I}$ do 7: Eligible drone types $\leftarrow \mathrm{list}$ of drone types that can serve $i$ 8: Cost for drone types in serving $i\leftarrow$ [] 9: for each drone type $d$ in eligible drone types do 10: Route index with minimum cost $\leftarrow -1$ 11: minimum cost $\leftarrow M$ 12: Route index with minimum remaining energy $\leftarrow -99$ 13: minimum remaining energy $\leftarrow M$ 14: ${\mathcal{L}}_{\mathrm{Ed}}^{\prime }$ , ${\mathcal{L}}_{\mathrm{Ed}}^{″}\leftarrow$ copy of ${\mathcal{L}}_{\mathrm{Ed}}$ 15: ${\mathcal{L}}_{\mathrm{rd}}^{\prime }\leftarrow$ copy of ${\mathcal{L}}_{\mathrm{rd}}$ 16: for each route in ${\mathcal{L}}_{\mathrm{rd}}$ do 17: $k\leftarrow$ last location in the current route 18: $\mathrm{batRepA}{t}_0$ , $\mathrm{batRepA}{t}_k$ , remaining energy $\leftarrow \mathrm{BatteryReplacement}$ ( $d,k,i$ , route index) 19: Update the remaining energy of the current route in ${\mathcal{L}}_{\mathrm{Ed}}^{\prime }$ 20: Potential cost $\leftarrow E_{k0}^d+E_{0i}^d+{\pi }_{\mathrm{bat}}^d\times (\mathrm{batRepA}{t}_0+\mathrm{batRepA}{t}_k)$ 21: Time window satisfied Boolean $\leftarrow \mathrm{CheckTimeWindow}$ ( $d,k,i$ , route index, $\mathrm{batRepA}{t}_0,\mathrm{batRepA}{t}_k$ ) 22: if Time window is satisfied and potential cost < minimum cost then 23: Minimum cost $\leftarrow \mathrm{Potential}$ cost 24: Route index with minimum cost $\leftarrow \mathrm{current}$ route index 25: if Remaining energy of the current route < minimum remaining energy then 26: Route index with minimum remaining energy $\leftarrow \mathrm{current}$ route index 27: if Route index with minimum cost $=-1$ then 28: Potential cost $\leftarrow {\pi }_d+E_{0i}^d$ 29: Add the potential cost to the list of costs for drone types 30: Add $[0,i]$ to ${\mathcal{L}}_{\mathrm{rd}}^{\prime }$ 31: else 32: Select the route, $r$ , with the minimum cost and minimum remaining energy 33: Add location $i$ to route $r$ in ${\mathcal{L}}_{\mathrm{rd}}^{\prime }$ 34: Update the remaining energy of route $r$ in ${\mathcal{L}}_{\mathrm{Ed}}^{″}$ 35: Add the cost of route $r$ to the list of costs for drone types 36: Selected drone type, $d^{*}\leftarrow$ select the drone type with the minimum cost from the list of costs for drone types 37: ObjectiveValue $\leftarrow \mathrm{ObjectiveValue}$ + cost for the selected drone type $d^{*}$ 38: ${\mathrm{nBat}}_d^{*}\leftarrow {\mathrm{nBat}}_d^{*}+\mathrm{batRepA}{t}_0+\mathrm{batRepA}{t}_k$ 39: Update ${\mathcal{L}}_{\mathrm{rd}}\leftarrow {\mathcal{L}}_{\mathrm{rd}}^{\prime }$ and ${\mathcal{L}}_{\mathrm{Ed}}\leftarrow {\mathcal{L}}_{\mathrm{Ed}}^{″}$ 40: return ${\mathcal{L}}_{\mathrm{rd}},{\mathrm{nBat}}_d^{*}$ , ObjectiveValue

Algorithm 2.

Battery Replacement Decisions and Time Window Check

Input $d,k,i$ , route index Output $\mathrm{batRepA}{t}_0,\mathrm{batRepA}{t}_k,\mathrm{rE}$ 1: procedureBatteryReplacement ( $d,k,i$ , route index) 2: $\mathrm{batRepA}{t}_k\leftarrow$ False 3: $\mathrm{batRepA}{t}_0\leftarrow$ False 4: $\mathrm{rE}\leftarrow$ get the remaining energy of current route from the list of remaining energies of drone type $d$ 5: if $\mathrm{rE}-E_k{0}^d<b_d^{\min }$ then 6: $\mathrm{batRepA}{t}_k\leftarrow$ True 7: $\mathrm{rE}\leftarrow b_d^{\max }-E_k{0}^d$ 8: else 9: $\mathrm{rE}\leftarrow \mathrm{rE}-E_k{0}^d$ 10: if $\mathrm{rE}-E_0{i}^d<b_d^{\min }$ then 11: $\mathrm{batRepA}{t}_0\leftarrow$ True 12: $\mathrm{rE}\leftarrow b_d^{\max }-E_0{i}^d$ 13: else 14: $\mathrm{rE}\leftarrow \mathrm{rE}-E_0{i}^d$ 15: return $\mathrm{batRepA}{t}_0,\mathrm{batRepA}{t}_k,\mathrm{rE}$ 16: procedureCheckTimeWindow( $d,k,i$ , route index, $\mathrm{batRepA}{t}_0,\mathrm{batRepA}{t}_k$ ) 17: $p_i^d=\max [e_i,p_k^d+t_{\ell }^d+t_u^d+t_{k0}^d+t_{0i}^d+t_{\mathrm{bat}}^d\times (\mathrm{batRepA}{t}_0+\mathrm{batRepA}{t}_k)]$ 18: if $p_i^d\le {\ell }_i$ then 19: return True 20: else 21: return False

Genetic Algorithm

This section presents the proposed customized GA to solve our MTMFDSP-RD. To reach a good quality solution faster, we applied Algorithm 1 to get a proper upper bound on the required number of drones for each drone type. The main challenge in implementing the GA is to produce feasible chromosomes in each generation, as we must assign an eligible drone type to each delivery/pickup location accounting for each drone type’s delivery range and package weight limitations. Moreover, each route (i.e., the sequence of delivery/pickup locations visited by a drone) of a drone of a particular type must be feasible relating to the release and due times of the delivery/pickup locations. As the fundamental nature of GA is very random, we must check each chromosome in every generation (i.e., the initial generation after generating the chromosomes, and the subsequent generations after applying crossover and mutation) and apply the necessary modifications to the infeasible genes to make feasible chromosomes. Procedure “ModifyChromosome” in Algorithm 3 shows the detailed steps of the necessary modifications to an infeasible chromosome. Algorithm 3 describes the detailed steps of the proposed GA.

As shown in Algorithm 3, after determining the sets of deliveries and the upper bound on the required number of drones for each drone type from the Algorithm 1, the chromosomes are generated for the initial population. Each chromosome contains $n$ genes as $n$ number of delivery/pickup locations. Genes are created randomly in a way that gene $i$ in each chromosome represents the drone type and drone number (i.e., index) used to serve the $i^{\mathrm{th}}$ delivery. For example, a gene can have the value of $A9$ representing that this delivery location is served with the ninth drone of drone type $A$ . The population size is set to be 10,000 to create the initial population. Then, the initial population is checked for feasibility, and necessary modifications are made based on Procedure “ModifyChromosome” in Algorithm 3. After conducting necessary modifications and making all the chromosomes feasible, the fitness function of each of the chromosomes is calculated, and the best fitness value is obtained among the entire population. To form the next generation of chromosomes, the selection method selects $b\mathrm{\%}$ of chromosomes with the best values to be directly transferred to the next generation. Then, $c\mathrm{\%}$ of the remaining chromosomes with the probability of $p$ is selected for uniform crossover to make new chromosomes for the next generation. Furthermore, the remaining chromosomes are mutated as a 1-point mutation, and the mutated chromosomes are also transferred to the next generation. As some of the chromosomes of the next generation are produced by crossover and mutation, it is possible for them to become infeasible due to the release and due times (i.e., time windows) and drone-type-related constraints. Therefore, the next generation is also checked for feasibility and is modified using the modification algorithm shown in Procedure “ModifyChromosome” of Algorithm 3. Then, we apply the same procedure for all of the generations until the algorithm terminates. To check the convergence of the algorithm, the average of the best values among the previous $k$ generations is calculated. If the best fitness function for the current generation is not improved by more than $a\mathrm{\%}$ compared to the average of the best values among the previous $k$ generations, the algorithm terminates and returns the best chromosome of the current generation as the best (i.e., final) solution. The parameters of the customized GA implemented in this paper are as follows: $a=0.1,b=1,c=90,k=4,p=0.5$ , which are selected by parameter tuning.

The modification procedure shown in Procedure “ModifyChromosome” of Algorithm 3 checks the input chromosome for the time window and drone-type-related constraints—drone range and package weight-carrying capacity. Firstly, each gene is checked to determine whether the assigned drone type can serve the associated delivery/pickup location or not based on the drone range and package weight-carrying capacity. If the assigned drone type cannot serve the associated location, a random drone number from the eligible drone types is randomly assigned to this gene. After checking all genes for the drone-type-related constraints, the delivery sequences in the chromosome are determined. The delivery/pickup locations that are served by the same drone number of a particular drone type are considered as a delivery sequence. This delivery sequence is checked for the release time and time window constraint satisfaction. If the delivery sequence is infeasible, the drone number of the first gene in this delivery sequence that violates the release time and time window constraint is changed to a new random drone number of the same drone type. This process is applied to the subsequent genes in this delivery sequence. This procedure is conducted for all the delivery sequences until all of them become feasible, and the resulting chromosome is the modified feasible chromosome.

Algorithm 3.

Genetic Algorithm

Input Delivery set $\mathcal{D}$ , drone types $\{A,B,C\}$ , earliest pick-up times, GA parameters Output Best Chromosome 1: procedureGeneticAlgorithm(Delivery set $\mathcal{D}$ , drone types $\{A,B,C\}$ , earliest pick-up times, GA parameters) 2: Sort deliveries in ascending order based on earliest pick-up time 3: Determine deliverability sets $S_1$ , $S_2$ , and $S_3$ for each drone type 4: Estimate upper bound on number of drones of each type using a Greedy Algorithm 5: Define GA parameters: population size, selection method, mutation rate, crossover rate, etc. 6: Define chromosome structure: each gene $\to \mathrm{drone}$ type and drone ID for corresponding delivery 7: Generate initial population 8: while termination criterion not met do 9: Call ModifyChromosome() to ensure feasibility 10: Evaluate fitness for all chromosomes 11: Compute best fitness value and average over past $k$ generations 12: if solution improved less than threshold $a\mathrm{\%}$ then 13: break 14: Select top $b\mathrm{\%}$ chromosomes for the next generation 15: Apply uniform crossover on $c\mathrm{\%}$ chromosomes with probability $p$ 16: Apply 1-point mutation on remaining chromosomes 17: Form new generation 18: Return Best Chromosome 19: procedureModifyChromosome(Chromosome) 20: for each gene $i$ in chromosome do 21: if delivery $i$ cannot be served by its current drone type then 22: Assign a random feasible drone type to gene $i$ 23: Determine all feasible delivery sequences for different drones 24: for each delivery sequence $t$ do 25: if sequence $t$ violates time window constraint then 26: Identify the first violating delivery in sequence $t$ 27: Replace its gene value with a randomly selected drone of the same type 28: Return Modified Chromosome

Numerical Results and Managerial Insights

In this section, we provide insights into the following research questions by applying our proposed mathematical model and the solution algorithms on the blood delivery logistics case study:

1. How do the runtime and solution quality of different algorithms compare to each other and how do model parameters affect their runtime?

2. How does using a mixed fleet of drones affect the total energy consumption and the total cost?

3. How do the required fleet composition and fleet size change with the pickup time window, package weight, minimum required energy in the battery, and the flight path?

4. How frequently must batteries of different drones be replaced as the package weight, flight path, and minimum required energy change?

5. How do the drone operating parameters—package weight, flight path, and the minimum required energy in the battery—affect the total energy consumption?

Case Study

This section presents a case study to evaluate the performance of the proposed MILP model and the algorithms in solving our MTMFDSP-RD for a blood sample delivery logistics problem. This case study, based on actual blood sample delivery and drone flight test data, is used to provide key managerial insights into using a mixed fleet of aerial drones for delivering time-sensitive medical items that would be beneficial for efficiently organizing and scheduling a fleet of drones in medical delivery logistics.

Blood Sample Delivery Data

The actual blood sample delivery data used in this case study are from the Interpath Laboratory, Inc. ( 46 )—a healthcare logistics company—located in Pendleton, Oregon, USA, which is a fully integrated, for-profit, clinical, and anatomic pathology medical laboratory. In Pendleton, the Interpath Laboratory provides services to 11 other healthcare business locations (e.g., hospitals/clinics) in the region (as shown in Figure 4a).

Figure 4.

Demonstration of delivery data and flight path: (a) the drone depot and the pickup locations, (b) a drone flight path with a circular no-fly zone. (Color online only.)

In this study, we considered the Interpath Laboratory as the central depot that hosts a mixed fleet of drones. Each drone flies from the depot (i.e., Interpath Laboratory location) to pick up the blood samples from other healthcare business locations and return to the depot. As the customers in this study are also other healthcare business locations (i.e., healthcare provider locations), they have battery swapping and charging infrastructure for drones. We used one week of customer orders that the Interpath Laboratory received from all 11 other healthcare business locations based on their needs. As the packages in this delivery system are blood samples, each order comes with an earliest pickup time (i.e., order release time) and a maximum delayed pickup time (i.e., due time) by which a drone must pick up the package from where it originates. This study considers the difference between these release and due times as the pickup time window that is different for different pickup tasks. The drone fleet operator knows each order’s location (i.e., latitude and longitude), timestamp, pickup time window, and package weight. A drone cannot reach the pickup location before the earliest pickup time, as the package is not ready before this timestamp. The original blood sample delivery data contain 81 pickup orders in a week of five business days. Based on these original data, we generated larger instances by increasing the frequency of the orders in each day to evaluate the performance of our MILP model and the proposed heuristic algorithms. It is evident from Figure 4a that the depot location is different from the geometric center of the healthcare business locations. Some pickup locations are close to the depot (blue marked locations in Figure 4a), whereas some other locations are far away (i.e., yellow marked locations in Figure 4a), as much as 30 mi (i.e., 48 km) from the depot. This real-world scenario necessitates using different types of drones to ensure an efficient delivery system as different drone types are suitable for delivery/pickup tasks at different distances due to their physical characteristics (e.g., speed, battery capacity, energy consumption). Besides, there exist different no-fly zones, such as airports, elementary schools, and municipal courts, in the flight path of drones. Some of these “no-fly zones” for the case study data are red-flag marked in Figure 4a that the drones need to avoid by detouring via some waypoints. Figure 4b demonstrates the possible flight paths (dashed line: straight path, blue line: road network, and red line: detouring due to no-fly zone) for a drone flying from a point A to another point B. (Color online only.)

Drone Types and Data

We observed that the healthcare business locations are at different distances from the depot and there are no-fly zones in the flight path of drones, which makes different drones suitable for pickup/delivery tasks at different healthcare business locations, necessitating the use of a mixed fleet of drones. Therefore, in this case study, we considered a mixed fleet of drones consisting of three drone types—quadcopter (Tarot 650), hexacopter (DJI Matrice 600 Pro), and fixed wing VTOL (Wingcopter 198)—for cost- and energy-efficient delivery of the blood samples. In their study, Bhuiyan et al. ( 13 ) used two types of rotary drones—quadcopter (Tarot 650), and hexacopter (DJI Matrice 600 Pro)—to investigate the drone-based delivery in a prepared food delivery case study. In our case study, along with these two rotary drones, we also considered VTOL (Wingcopter 198) as the rotary drones cannot serve some of the healthcare business locations (yellow marked locations in Figure 4a) due to delivery range limitations and package weight carrying capacity restrictions. VTOL has a longer delivery range and package weight carrying capacity than rotary drones and can serve every healthcare business location in the blood sample logistics case study. However, VTOL has a higher energy consumption rate in ascending and descending while having a lower energy consumption rate in forward flight for longer distances compared to rotary drones. Thus, it makes this drone type less energy efficient in delivering/picking up packages closer to the depot, but more efficient than the rotary drones in delivering/picking up packages at longer distances.

In the numerical experiments, we used the drone power consumption and flight time data collected from our drone flight tests conducted at the Idaho National Laboratory (INL) test site. We obtained the power consumption and flight time of different drone types from the flight test data for different values of the parameters—drone speed, package weight, and flight path (i.e., straight versus with turn). We used an ML method (discussed in the section “Estimation of energy consumption of drones”) on these experimental data to estimate the power and energy consumption of drone types for a range of parameter values (e.g., different values of package weight) used in our blood delivery case study. The battery capacity, speed, package weight carrying capacity, and costs of the three drone types in all numerical experiments of this paper are taken from Bhuiyan et al. ( 13 ).

Experimental Setup

To solve the proposed MILP model in the section “MILP model”, we implemented the model using Python 3.10 with CPLEX optimizer 22.10 ( 47 ). We implemented our proposed GA and GH algorithm using Python 3.10. To test the computational performance of the model and the algorithms, we synthetically generated larger problem instances based on the original blood sample delivery data by increasing the frequency of customers’ orders in each day over one week of planning horizon. We considered three different package weight distributions in the numerical experiments to gain insights into the effect of package weights on fleet size, fleet composition, and total energy consumption. We refer to these distributions as distribution A, distribution B, and distribution C. Distribution A is the original distribution of package weights in the blood sample logistics data from the Interpath Laboratory, Inc., a healthcare logistics company, where 88% of the package weights lie between 0.5 to 1 kg, 6% between 1 to 1.5 kg, and 6% between 1.5 to 2 kg. We used this package weight distribution for all the sensitivity analysis. We synthetically generated distributions B and C. Distribution B includes 33.33% of the packages weighing between 0.5 to 1 kg, 33.33% between 1 to 1.5 kg, and 33.33% between 1.5 to 2 kg. Distribution C includes 15% of the packages weighing between 0.5 to 1 kg, 15% between 1 to 1.5 kg, 10% between 1.5 to 2 kg, 20% between 2 to 2.5 kg, 20% between 2.5 to 3.5 kg, and 20% between 3.5 to 5.5 kg. These distributions and the corresponding upper and lower bounds of the package weights are determined based on our discussion with the Interpath Laboratory about their historical blood sample logistics data. We also conducted numerical experiments for different combinations of the parameters shown in Table 2 to analyze the effects of these parameters on total energy consumption, the required fleet size, required fleet composition, and the required number of battery replacements.

Table 2

Parameters for Sensitivity Analysis

Parameters	Values
Minimum required battery energy, $b_d^{\min }$ (% of initial battery energy)	10, 15, 20, 25
Pickup time window (min)	30, 45, 60, 120, 180

Runtime of the Solution Methods

In this section, we present the computational performance of the proposed MILP model, valid inequality, GH algorithm, and the GA. Table 3 presents the (1) runtime of the CPLEX solver in solving the MILP model with and without the valid inequality, (2) runtime and solution quality (optimality gap) of the GH algorithm, and (3) runtime and solution quality (optimality gap) of the GA. The “Basic”, “Accelerated”, “GH”, and “GA” columns in Table 3 denote the CPLEX solver without the valid inequality, the CPLEX solver with the valid inequality, the GH algorithm, and the GA, respectively. We conducted computational experiments by varying the size of the blood delivery logistics network (i.e., delivery data) and the pickup time window at the healthcare business locations. The smallest network size (i.e., 81 deliveries) is the original blood delivery logistics data obtained from Interpath laboratory, Inc., a healthcare logistics company, whereas the largest network size (i.e., 405 deliveries) is 5 times the original data from Interpath laboratory.

Table 3.

Runtime of the Solution Procedures

Delivery data (number of deliveries)	Time window (min)	Runtime (s)				Optimality gap (%)
		Basic	Accelerated	GH	GA	GH	GA
81	30	*	1161	0.01	1087	5.02	1.1
	45	*	1279	0.01	777	5.09	0.12
	60	*	1533	0.02	698	11.5	13.8
	120	*	10251	0.01	438	3.61	2.47
	180	*	*	0.01	755	na	na
243	30	*	*	0.03	5223	na	na
	45	*	*	0.04	3363	na	na
	60	*	*	0.05	2970	na	na
	120	*	*	0.04	1524	na	na
	180	*	*	0.05	1844	na	na
405	30	*	*	0.07	*	na	na
	45	*	*	0.08	*	na	na
	60	*	*	0.09	*	na	na
	120	*	*	0.08	*	na	na
	180	*	*	0.09	*	na	na

Note: GA = genetic algorithm; GH = greedy heuristic algorithm. * = parameter combinations that cannot be solved within the three-hour time limit; na = not applicable (as the exact solution is not obtained, the corresponding optimality gap cannot be computed).

Table 3 demonstrates that CPLEX cannot solve the MILP model without the valid inequality in three hours. The parameter combinations that cannot be solved within the three-hour time limit are marked with an asterisk (*) in Table 3. The valid inequality improves the computational efficiency and enables CPLEX to solve the model for the original data (i.e., 81 pickup operations) for most of the pickup time windows. The runtime of CPLEX increases as the pickup time window becomes more relaxed. This is because a wider time window allows more possible routes for the drones and thus the optimization routine needs to evaluate a larger number of candidate solutions, resulting in an increased runtime. However, CPLEX (including the valid inequality) cannot solve the model for delivery data larger than the original one. We see that the network size (i.e., size of delivery data) has a negligible effect on the runtime of the GH algorithm, whereas the runtime of GA increases as the network size increases. As the network size increases, the chromosome size also increases in GA. This increases the frequency of chromosome modification (to make infeasible chromosomes feasible) significantly, resulting in a longer time for the GA to converge. However, the GH algorithm does not suffer from this limitation as GA does. The runtime of the GA has a decreasing trend, despite having some randomness, as the pickup time window increases. This is because, as the time window becomes stricter, more routes (chromosomes) become infeasible necessitating more chromosome modifications and thus more iterations, eventually increasing the GA runtime. However, as the time window increases, less routes (chromosomes) are infeasible, reducing the need for chromosome modifications and thus the runtime. The solution quality of the GA is slightly better than the GH algorithm in most of the problem instances.

Effect of Mixed Fleet of Drones on Total Costs and Energy Consumption

Among the three different types of drones used in this study, the rotary quadcopter drone (Tarot 650) has the least delivery range and package weight carrying capacity, whereas the VTOL (Wingcopter 198) has the largest package weight carrying capacity and longest delivery range. We observed from our drone flight tests that the rotary quadcopter drone is the most energy efficient among the three drone types in delivering a package to the same limited distance, whereas the VTOL is more energy efficient in delivering packages to a longer distance. As the healthcare business locations are at different distances from the depot in our case study that also reflects the real-world delivery scenario, different drone types are suitable for delivery/pickup tasks at different healthcare business locations, reflecting the potential benefit of using a mixed fleet of drones in the logistics problem. Therefore, we evaluated the effect of using a mixed fleet of three drone types on the total energy consumption and total cost of the fleet in delivering blood samples.

Figure 5 demonstrates the variation in the total cost as the available number of drone types in the mixed fleet varies. We see that using a mixed fleet consisting of all three types of drones results in a lower total cost compared to using two drone types (i.e., hexacopter and VTOL) and one drone type (i.e., a homogeneous drone fleet with only VTOL drones). This is because the rotary quadcopter (i.e., Tarot 650) is the cheapest drone type for shorter-range deliveries, whereas the VTOL (i.e., Wingcopter 198) is the most expensive drone type. Therefore, the mixed fleet optimization uses the quadcopter and hexacopter to conduct pickup tasks from the healthcare business locations within their range, whereas it uses VTOL for the healthcare business locations outside the range of rotary drones or when the package weight exceeds the capacity of rotary drones. The numerical results also demonstrate that the total energy consumption decreases as the available number of drone types in the mixed fleet increases. Specifically, using a mixed fleet comprising rotary quadcopter, rotary hexacopter, and VTOL decreases the total energy consumption by 28.7% compared to a homogeneous fleet of only VTOL drones. This is because VTOL consumes much higher energy in ascending and descending compared to the rotary drones, which makes VTOL the most energy-intensive drone type for shorter-range deliveries. However, VTOL drones consume less energy than the rotary drones while flying straight, making VTOL drones more suitable for longer-distance deliveries compared to rotary drones.

Figure 5.

Effect of mixed fleet of drones on total costs.

It is evident from these results that the business owners can substantially reduce the total costs and total energy consumption by using a mixed fleet of drones consisting of these three types. For instance, for a pickup time window of 180 min, using a mixed fleet consisting of all three types of drones reduces the total costs and total energy consumption by 18.18% and 28.7%, respectively, compared to using only VTOL drones in the fleet.

Effect of Package Weight Distribution and Drone Flight Path on Fleet Size and Fleet Composition, Required Number of Additional Batteries, and Energy Consumption

We conducted numerical experiments to investigate the effect of package weight distribution and drone flight path (i.e., straight path versus no-fly zone detouring and over the road network) on the required fleet size and fleet composition (i.e., number of drones of each type in the fleet), required number of additional batteries to operate the fleet, and total energy consumption. In this analysis, we used the different package weight distributions (discussed in the section “Experimental setup”) generated from the original blood sample delivery logistics data. We used the ML model (presented in the section “Estimation of energy consumption of drones”) and the drone flight test data to estimate the energy consumption of different drone types for these package weight distributions. We used the Open Street Routing Map (OSRM) Application Programming Interface (API) to compute the road distance and the number of turns a drone takes while flying to each healthcare location from the depot and then used these distances, turns, and the estimated energy consumption using ML in the proposed MILP model as input parameters. In the numerical experiments in this subsection, we used a 30-min pickup time window for the healthcare business locations. Additionally, we set the minimum required battery energy for each drone type at 15% of its initial battery capacity. Table 4 demonstrates the effect of package weight distribution and flight path on the required fleet size and fleet composition, total energy consumption, and the required number of additional batteries.

Table 4.

Effect of Package Weight Distribution and Drone Flight Path on Fleet Size and Fleet Composition, Total Energy Consumption, and the Required Number of Additional Batteries

Package weight distribution and flight path	Required number of drones			Total energy consumption (MJ)	Required number of batteries
	Quad	Hex	VTOL		Quad	Hex	VTOL
Distribution A
Straight	1	0	4	98.3	10	0	20
NFZ	1	0	4	104.6	13	0	21
RN	1	1	4	131.3	2	4	31
Distribution B
Straight	0	1	4	106.2	0	10	18
NFZ	0	1	4	112.8	0	8	24
RN	1	1	4	138.1	2	5	32
Distribution C
Straight	0	1	4	114.7	0	4	23
NFZ	1	1	4	119.2	2	5	25
RN	1	2	4	140.5	5	4	28

Note: VTOL = vertical-takeoff-and-landing; NFZ = no-fly zone for flight path; RD = road network for flight path; Quad = rotary quadcopter; Hex = rotary hexacopter.

We see from Table 4 that the mixed fleet includes drones with larger package weight carrying capacity and higher battery capacity (i.e., using hexacopter and VTOL than quadcopter) as the package weight distribution includes heavier packages. This is because, in the mixed fleet of drones, hexacopter and VTOL drone types have larger package weight carrying capacity and battery capacity than quadcopter. Therefore, the required fleet composition includes more hexacopter and VTOL drones as the package weight distribution of the healthcare business locations includes a larger number of heavier packages.

We also see from Table 4 that the required fleet size increases and fleet composition changes as the drones need to detour through waypoints and fly over the road networks to avoid no-fly zones compared to flying in a straight path. This is because, with detours and flying over the road networks, a drone needs to fly for a longer distance to travel between two points in the network compared to flying in a straight path. Due to this additional time while detouring and flying over the road networks, it is not possible to maintain the same pickup time window with the same number of drones as it is possible while flying in a straight path, resulting in an increased number of drones in the fleet. For a 30-min pickup time window, on average over the three different package weight distributions, the required fleet size increases by 6.67% and 26.67% as the drones in the mixed fleet detour through waypoints and fly over the road networks, respectively, compared to flying in a straight path.

Moreover, we see from Table 4 that the total energy consumption by the mixed fleet of drones increases in each package weight distribution as the drones detour through waypoints and fly over the road networks to avoid no-fly zones compared to flying in a straight path. This is because while detouring through waypoints and flying over the road networks, drones need to fly for a longer distance to travel between two points in the network compared to flying in a straight path, resulting in an increased energy consumption. Notably, the total energy consumption, on average, increases by 28.25% and 5.31% when drones fly over the road networks and detour through waypoints, respectively, compared to flying in a straight path. Furthermore, we see from Table 4 that the total energy consumption increases as the package weight distribution includes more heavier packages, which is intuitive as the package weight significantly affects the drone energy consumption.

Additionally, we see from Table 4 that the required number of battery replacements increases as the drones need to fly with detours through waypoints and over the road networks to avoid no-fly zones compared to flying in a straight path for each package weight distribution. This is because, flying with detours through waypoints and over the road networks results in a significantly higher energy consumption by the drones than flying in a straight path. This increased energy consumption in each pickup/delivery task results in reaching the minimum required battery energy more frequently for the drones, and thus the frequency of battery replacements increases. This analysis reveals that, on average over three different package weight distributions, battery replacement increases by 15.32% and 32.96% when drones detour through some waypoints and traverse over road networks, respectively, instead of flying in a straight path.

However, we see from Table 4 that the total number of battery replacements does not follow any particular increasing/decreasing trend as the package weight distributions vary. This is because, in the package weight distributions with more heavier packages, the optimal fleet composition consists of more larger drones (i.e., hexacopter and VTOL than the quadcopter) that have larger battery capacity and thus require fewer battery replacements. For instance, in the case of “straight” flight path, moving from package weight distribution A to B, all the quadcopters in distribution A are replaced by hexacopters in distribution B, eventually resulting in a reduced number of total battery replacements as the hexacopter’s battery capacity is more than three times the battery capacity of quadcopter.

Effect of Time Window on Fleet Size and Fleet Composition

We conducted numerical experiments to evaluate the effect of pickup time window requirements of the healthcare business locations on the required size and composition of the mixed fleet. Table 5 shows that the required fleet size decreases as the time window becomes more relaxed. As mentioned earlier, our problem is a multi-trip mixed-fleet drone scheduling problem with release and due times (MTMFDSP-RD), where we aim to have each drone perform multiple trips between the depot and the healthcare business locations over the planning horizon to reduce the required number of drones in the fleet in performing all the pickup tasks maintaining the pickup time window. Therefore, with a more relaxed pickup time window requirement, each drone can perform more pickup tasks from the healthcare business locations, reducing the required number of drones in the fleet to perform all the pickup tasks over the planning horizon. For package weight distribution A, the required fleet size (i.e., required number of drones) increases by 100% as the pickup time window decreases from 180 min to 30 min.

Table 5.

Time Window on Fleet Size and Fleet Composition

Time window (min)	Required number of drones
	Quadcopter	Hexacopter	VTOL
30	1	1	4
45	1	1	3
60	1	0	3
120	1	0	2
180	1	0	2

Note: VTOL = vertical-takeoff-and-landing.

Table 5 also demonstrates how the required fleet composition changes as the required pickup time window varies. We see that the strict time windows require using all three types of drones, whereas the smallest and the largest drones are used as the time window becomes more relaxed. For instance, the 30 and 45 min pickup time window requires using all three types of drones. The reason for this phenomenon is the model’s tendency to prioritize using shorter-range drones (i.e., quadcopters) due to their high energy efficiency and low cost, and then use the longest-range drones (i.e., VTOL) to perform pickup tasks outside the range of quadcopters. However, the battery capacity of these quadcopters is very limited and therefore needs to be replaced more frequently compared to the larger drones (i.e., hexacopter and VTOL). As each battery swap (replacement) takes several minutes, more frequent battery replacements reduce the drone utilization and require a larger number of drones to satisfy the pickup time window while using quadcopters. Therefore, the proposed optimization model decides to use some hexacopters having much larger battery capacity than quadcopters instead of using many quadcopters to satisfy the strict pickup time window requirement. However, when the pickup time window is more relaxed, such as 60 to 180 min, despite more frequent battery replacements, quadcopters still can meet the time window requirement. However, the drone fleet still requires using VTOL to perform pickup tasks for healthcare business locations outside the delivery range of quadcopters. Additionally, with a strict time window (e.g., 30 min), a larger number of VTOL drones are used in the mixed fleet, due to their higher speed, to maintain the required time window compared to the case of a more relaxed time window.

Effect of Minimum Required Battery Energy on the Required Number of Additional Batteries

As previously discussed, we must maintain a minimum required energy in the drone battery upon arrival at any node in the delivery network to ensure the safe operation of the drones. If the remaining energy in a drone’s battery falls below the minimum threshold value, it can fall down during flight. In most of our experiments, we considered this minimum required battery energy to be 15% of the initial battery energy. However, practitioners may use a higher threshold for the minimum energy requirement to be more risk-averse in operating their drone fleet. We find that the minimum required battery energy significantly affects the frequency of the required number of battery replacements (i.e., the required number of additional batteries) in the pickup operations from the healthcare business locations. We conducted numerical experiments for four distinct values of the minimum required battery energy to evaluate the effect of the minimum required battery energy on the required number of additional batteries, as shown in Table 6. In these experiments, we used the package weight distribution of the original blood delivery logistics data (i.e., distribution A) and a 30-min pickup time window.

Table 6.

Effect of Minimum Required Battery Energy on the Required Number of Additional Batteries

Minimum required battery energy (% of initial battery energy)	Required number of additional batteries
	Quadcopter	Hexacopter	VTOL
10	9	0	17
15	10	0	20
20	14	5	21
25	8	10	26

Note: VTOL = vertical-takeoff-and-landing.

Table 6 demonstrates that a more conservative minimum required battery energy forces the model to use larger drones (i.e., hexacopter and VTOL), and thus requires more frequent battery replacements for larger drones compared to the smaller drones. This phenomenon is intuitive because a conservative minimum required battery energy (i.e., 25% of the initial battery energy) makes more healthcare business locations out of range for smaller drones (i.e., quadcopters) which could be reached by these smaller drones if the minimum required energy is set to 10% or 15%. Thus, the model must assign larger drones (i.e., hexacopters or VTOLs) to serve these healthcare business locations, resulting in more frequent battery replacements for these larger drones. Increasing the minimum battery energy requirement results in reaching the minimum threshold more frequently during drone operations, thereby increasing the required number of additional batteries (i.e., frequency of battery replacements). Table 6 shows a significant rise (i.e., 69.23%) in the required number of battery replacements as the minimum required energy increases from 10% to 25% of the initial battery energy.

Comparison with Homogeneous Fleet and Pre-Determined Battery Replacement

We conducted numerical experiments to compare the benefit of our proposed MTMFDSP-RD in leveraging a mixed fleet of drones and an endogenous battery replacement strategy with the existing literature that uses a homogeneous fleet of drones and pre-determined battery replacement at each node. In this analysis, we used three different drone types, including rotary quadcopter (Tarot 650), rotary hexacopter (DJI Matrice 600 Pro), and VTOL (Wingcopter 198) in the mixed fleet. In the homogeneous fleet, we considered VTOL (Wingcopter 198) that can deliver packages to all locations. Table 7 shows the number of healthcare business locations in our case study that are out of the delivery range of each drone type. Moreover, Table 8 demonstrates the effect of available drone types in a fleet (i.e., a homogeneous drone fleet with only VTOL drones versus a mixed fleet comprising rotary quadcopter, hexacopter, and VTOL) on the total cost and energy consumption.

Table 7.

Number of Healthcare Business Locations out of the Delivery Range of Different Drone Types

Drone type	Healthcare business locations out of the delivery range
	Number	Percentage of total healthcare business locations (%)
Quadcopter	22	27
Hexacopter	19	23
VTOL	0	0

Note: VTOL = vertical-takeoff-and-landing.

Table 8.

Effect of Mixed Fleet of Drones on Total Costs and Energy Consumption

Fleet type	Time window (min)	Total cost ( $\times \mathrm{\$}1,000$ )	Total energy consumption (MJ)
Homogeneous (only VTOL)	30	335	135.0
	45	274	131.0
	60	214	127.4
Mixed fleet (quadcopter, hexacopter, and VTOL)	30	296	96.5
	45	244	93.2
	60	178	90.7

Note: VTOL = vertical-takeoff-and-landing.

Table 7 demonstrates that using only rotary drones—quadcopter (Tarot 650) and hexacopter (DJI Matrice 600 Pro)—in the fleet results in a considerably larger number of healthcare business locations out of the delivery range of drones. Specifically, a homogenous fleet comprising rotary quadcopters results in 27% of healthcare business locations out of the delivery range, while using only rotary hexacopters in the drone fleet makes 23% of healthcare business locations out of the delivery range of drones. In contrast, VTOL (i.e., Wingcopter 198 in our case study) is the only drone type in our drone fleet that can serve all healthcare business locations.

However, Table 8 shows that using only VTOL in the fleet increases the total cost and energy consumption by 15% and 40%, respectively, on average across different time windows compared to a mixed fleet of drones. This is because, VTOL (Wingcopter 198) has a significantly higher cost than the rotary drones—quadcopter (Tarot 650) and hexacopter (DJI Matrice 600 Pro). However, VTOL has faster speed, longer delivery range, and larger package weight carrying capacity than rotary drones. Also, VTOL consumes much higher power while ascending and descending, whereas it consumes lower power during forward flight compared to rotary drones, making VTOL suitable for longer-distance delivery with heavier packages. On the other hand, the rotary drones are limited to shorter-range delivery due to smaller battery capacity, slower speed, and less package weight carrying capacity, and thus are unable to deliver packages to far-away locations. Therefore, VTOL is more cost- and energy-efficient in delivering heavier packages at longer distances from the depot, whereas the rotary drones are more efficient in delivering lighter packages closer to the depot. Using a mixed fleet allows leveraging the relative benefits of different drone types to make the logistics operation more cost- and energy-efficient than a homogeneous fleet.

We also conducted numerical experiments to compare our proposed endogenous battery replacement strategy with the pre-determined battery replacement strategy at each node in the transportation network, as shown in Table 9. In this analysis, we used the mixed drone fleet comprising rotary quadcopter, rotary hexacopter, and VTOL. In the pre-determined battery replacement strategy, we replace the drones’ batteries at each node (i.e., both depot and healthcare business locations) in the delivery network.

Table 9.

Effect of Battery Replacement Strategy on the Total Cost, Required Number of Drones, and Required Number of Additional Batteries

TW (min)	Pre-determined battery replacement				Endogenous battery replacement
	Cost ( $\times \mathrm{\$}1,000$ )	No. of drones	ADPD	No. of batteries	Cost ( $\times \mathrm{\$}1,000$ )	No. of drones	ADPD	No. of batteries
30	333	7	11.6	155	296	6	13.5	29
45	273	6	13.5	156	244	5	16.2	32
60	213	5	16.2	157	178	4	20.2	33

Note: TW = time window; ADPD = average number of healthcare business locations served by each drone.

Table 9 demonstrates that, compared to our proposed endogenous battery replacement strategy, on average across all time windows, replacing the battery of the drones at each node of the transportation network requires a 399% larger number of additional batteries to support the drone-based delivery operation. In other words, our proposed endogenous battery replacement strategy significantly reduces the required number of additional batteries compared to the existing pre-determined battery replacement strategy. Additionally, unnecessarily replacing the batteries of drones at each node in the network increases the idle time of the drones, increasing the delay between consecutive delivery/pickup tasks. Therefore, compared to replacing the battery of drones in an as-needed fashion, pre-determined battery replacement at each node reduces the number of packages each drone can deliver while satisfying the release and due time of medical items, decreasing the utilization of drones. Table 9 shows that the average number of healthcare business locations served by each drone decreases by 83.1%, on average across all time windows, when replacing the battery of drones at each node compared to our proposed endogenous battery replacement strategy. Therefore, compared to the endogenous battery replacement strategy, replacing the battery at each node results in a 14.7% higher total cost and 20.6% higher required fleet size, which reflects the substantial cost-benefit offered by the proposed endogenous battery replacement strategy.

Practical Implications

In this paper, we proposed a new MILP model and faster algorithms for leveraging a mixed fleet of drones with distinct characteristics (e.g., cost, battery capacity, power consumption, package weight carrying capacity, and speed) and endogenous drone battery replacement to deliver time-sensitive medical items having distinct release and due times. We evaluated our model and solution algorithms based on a blood sample logistics case study from Interpath Laboratory, Inc.—a healthcare logistics company—and actual drone data using three drone types—rotary quadcopter (Tarot 650), rotary hexacopter (DJI Matrice 600 Pro), and VTOL (Wingcopter 198). Our proposed model and solution algorithms are general and can be applied to make up as well as efficiently route and schedule mixed drone fleets for other time-sensitive items delivery, such as whole blood, defibrillators, vaccines, and prepared food delivery. Practitioners and relevant companies can use our model and algorithms with their specific drone fleet and delivery data. Our solution algorithms can serve as fast and accurate decision-support tools for practitioners.

Our numerical results based on actual blood sample delivery data and drone data provide key managerial insights for healthcare business owners and medical logistics companies. Our numerical results highlight the benefit of using a mixed fleet of drones compared to a homogeneous drone fleet. Our proposed model and solution algorithms equip practitioners to evaluate the benefit of a mixed fleet of drones over a homogeneous fleet, and the benefit of endogenous battery replacement based on their specific drone and delivery data. We also studied the effect of different parameters, including pickup time window, package weight distribution, minimum required energy in the drone battery, and flight path on the required drone fleet composition and fleet size, required number of additional batteries, and the total energy consumption of the drone delivery system. These numerical results provide valuable insights for practitioners to help them understand how much the optimal routing and scheduling solutions change by varying each key parameter of the drone-based healthcare delivery system.

Moreover, as the main goal of the proposed MTMFDSP-RD is to minimize the cost of the required fleet size and the required number of additional batteries, our proposed model and solution algorithms are particularly beneficial for developing countries, where the budget is limited. Specifically, in cases of the high cost of specialized drone batteries, our proposed model with endogenous battery replacement is well-suited for developing countries, where budget constraints are significant and logistics system owners aim to find a delivery network with a minimum cost for collecting/delivering time-sensitive medical items from/to geographically distributed healthcare business locations.

Conclusion

We studied a multi-trip mixed-fleet aerial drone scheduling problem with release and due times, and endogenous battery replacements (MTMFDSP-RD)—timing and location of drone battery replacement is decided in the optimization routine accounting for the routing, remaining battery energy, and the minimum required battery energy—for time-sensitive medical supplies delivery. This problem seeks to optimally schedule and route a drone fleet containing different types of drones to determine the optimal fleet composition, required fleet size, required number of additional batteries, and total energy consumption. We used a supervised ML model to estimate the energy consumption of different types of drones for different operating parameters using actual drone flight test data. Introducing a new “route” concept, we proposed a new MILP to efficiently model the MTMFDSP-RD and introduced a valid inequality to improve the computational efficiency of the MILP model. We presented a new GH algorithm and a customized GA to solve the MTMFDSP-RD for large problem instances faster.

We evaluated the performance of the proposed MILP model and the solution algorithms in mixed fleet drone scheduling and routing with release and due times using an actual blood sample delivery case study from a healthcare logistics company (i.e., Interpath Laboratory, Inc.) located in Pendleton, Oregon, USA. We conducted numerical experiments using the actual drone flight test data and the actual blood sample delivery data to provide insights into (1) the computational performance of the proposed MILP model, valid inequality, GH algorithm, and the GA in solving MTMFDSP-RD for medical supplies delivery; and (2) the effect of drone fleet type and the drone operating parameters on the fleet composition, fleet size, additional batteries, and total energy consumption. Computational experiments demonstrate that the GH algorithm substantially outperforms the accelerated CPLEX (i.e., CPLEX with the valid inequality) and the GA while sacrificing the solution quality by a small amount. The solution quality of the GA is slightly better than the GH algorithm for some problem instances. While the network size (i.e., delivery data) has a negligible effect on the runtime of the GH algorithm, the runtime of GA increases significantly as the network size increases; the GA is not able to solve the largest network, which is 5 times the original blood delivery data in the case study.

Results demonstrate that the total cost and total energy consumption decrease as the mixed fleet contains all the drone types (i.e., quadcopter, hexacopter, and VTOL) compared to either using two drone types (i.e., hexacopter, and VTOL) or only VTOL drones; using a mixed fleet of three drone types reduces the total cost and total energy consumption up to 18.18% and 28.7%, respectively, compared to using a homogeneous fleet of only VTOL drones. The required fleet size, required number of additional batteries, and total energy consumption increase as the package weight distributions contain more heavier packages and drones fly with detouring through waypoints and over the road networks to avoid no-fly zones. For a 30-min pickup time window, on average, the required fleet size increases by 6.67% and 26.67% as the drones in the mixed fleet detour through waypoints and fly over the road networks, respectively, compared to flying in a straight path. Flying with detouring through waypoints requires traveling longer distance and time, which makes it impossible to maintain the same time window with the same number of drones as it is possible while flying straight. The required number of additional batteries to operate the mixed fleet increases, on average, by 15.32% and 32.96% as the drones detour through waypoints and fly over the road networks, respectively, compared to flying in a straight path.

The optimal composition of the mixed fleet includes larger number of drones with higher package weight-carrying capacity and battery capacity as the package weight distribution includes a larger number of heavier packages and drones’ flight path includes detouring. Results demonstrate that the required fleet size increases as the pickup time window of the healthcare business locations becomes stricter. For the original blood sample logistics data, the required fleet size increases by 100% as the pickup time window decreases from 180 min to 30 min. The required mixed fleet composition changes as the time window varies. Results demonstrate that all the three types of drones are used in the mixed fleet with strict time windows (i.e., 30 and 45 min) to satisfy the time window of the pickup tasks at the healthcare business locations. However, with a more relaxed time window, the mixed fleet mostly contains the smallest drones (also the most cost- and energy-efficient) and a few longest-range drones (i.e., VTOL) to pick up packages from healthcare business locations outside the range of smaller drones. Additionally, the mixed fleet contains a relatively larger number of drones with the highest speed when the time window is strict compared to a more relaxed time window.

The practical implication of this study is that practitioners can use our proposed optimization model and algorithms to schedule and route their mixed fleet of drones to determine the optimal fleet composition, required fleet size, required number of additional drone batteries, and energy efficiency of the fleet operation for their specific time-sensitive medical supplies delivery needs. A potential future extension of this paper is to study the mixed fleet of drones scheduling under uncertainties in demand and weather. Specifically, weather (i.e., wind, temperature) significantly affects the drone energy consumption and operability as some drone types cannot be used in adverse weather conditions. This problem can be modeled as a stochastic program to decide robust fleet composition and fleet size against uncertainty.