University of Southern California Dissertations and Theses

Computational geometric partitioning for vehicle routing

(USC Thesis Other)

Computational geometric partitioning for vehicle routing

PDF

Download a page range

Download transcript

Copy asset link

Request this asset

Transcript (if available)

Content COMPUTATIONAL GEOMETRIC PARTITIONING FOR VEHICLE ROUTING by Han Yu A Dissertation Presented to the FACULTY OF THE USC GRADUATE SCHOOL UNIVERSITY OF SOUTHERN CALIFORNIA In Partial Fulfillment of the Requirements for the Degree DOCTOR OF PHILOSOPHY (INDUSTRIAL & SYSTEM ENGINEERING) December 2023 Copyright 2023 Han Yu Dedication To everyone. ii Acknowledgements I would like to express my deep gratitude to all those who have supported me throughout this PhD journey, enabling me to reach its successful completion. I am profoundly grateful to my supervisor, John Carlsson, and Maged Dessouky, for their insightful guidance and unwavering support. I am thankful for their patience, encouragement, and the countless hours spent reviewing and refining my work. Also thanks to Sheng Liu, Stanley Lim, Sze-chuan Suen, and Fernando Ordonez for their invaluable guidance in my research projects. I would also like to convey my gratitude to the distinguished members of my dissertation and qualifica- tion committee, Petros A Ioannou and Andres Gomez. Their constructive feedback, valuable suggestions, and thought-provoking discussions have been crucial in refining the direction and quality of this research. I am incredibly grateful for my family and friends, whose support has meant the world to me. I want to express my deep appreciation to my friends at USC – Fangning Zheng, Yijing Lu, Ziyu He, Yanning Zhou, and Irene Li – for their unwavering encouragement and companionship. I also want to thank my friends outside of USC, Luying Sun, Kun He, and Shutian Li, for being there for me. I can’t forget to mention my companions in skiing and climbing – Emma Zhao, Tianbo Xie, Xiaorui Yang, Dailu Chen, Ranran, Ximeng, and Tianchen Yuan. Your friendship has been a continuous source of positivity in my life. To all my friends in China, your support during my PhD journey has been absolutely invaluable. I would like to acknowledge the assistance and insights of the research group, Suyanpeng Zhang, Yuanxiang Wang, Weizhi Lin, Yue Yu, Peng Dai, Bo Jo, their thoughtful conversations and collaboration iii enriched the research process. Their contributions played a crucial role in shaping the concepts and ideas presented in this dissertation. Lastly, I would like to express my appreciation to Xintao Zhao. His patience and support helped a lot in my PhD. In conclusion, this dissertation would not have been possible without the collective support of the individuals mentioned above. Their contributions have left an indelible mark on my academic and personal growth, for which I am truly grateful. iv TableofContents Dedication ii Acknowledgements iii List of Tables vii List of Figures viii Abstract x Chapter 1: Introduction 1 1.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1 1.2 Related Literature . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2 1.2.1 Region Partitioning in Logistics . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3 1.2.2 Stochastic and Dynamic Vehicle Routing . . . . . . . . . . . . . . . . . . . . . . . . 4 1.2.3 On-Time Delivery Operations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4 1.3 Application of Geometric Districting . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 Chapter 2: Convex Single-Depot Region Partitioning 7 2.1 Provably Good Region Partitioning . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8 2.1.1 Heterogeneous Service Times . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16 2.2 Optimal Partitioning Algorithms . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18 2.2.0.1 DividingR with a line . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21 2.2.0.2 Finding a feasible triple . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22 2.2.0.3 Recursive partitioning . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25 2.3 Numerical Experiments . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26 2.3.1 Benchmark Policies . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26 2.3.2 Simulation Study . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27 2.3.3 Uniform density function. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 28 2.3.4 Nonuniform density function. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30 2.3.5 A Real-World Case Study . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 33 Chapter 3: Orthogonal Partition 37 3.1 Optimal Partition . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 38 3.2 Algorithm . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 43 3.3 Numerical Analysis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 46 3.3.1 Experiment result at Los Angeles County . . . . . . . . . . . . . . . . . . . . . . . 47 3.3.2 Experiment result at Chicago . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51 v Chapter 4: Conclusions 55 Bibliography 56 vi ListofTables 2.1 Average delivery time (in minutes) of different fleet configurations and region partitioning policies in the case study. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 36 3.1 Average TSP distance and standard deviation (std) among subregions for convex and orthogonal partitions at Los Angeles. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 50 3.2 Average TSP distance and standard deviation (std) among subregions for convex and orthogonal partitions at Chicago. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 52 vii ListofFigures 2.1 A recursive partition whereR is a square, µ 1 is the density indicated by the shaded contour map, andµ 2 is uniform; the end result in 2.1e partitionsR intom = 13 pieces, all of which have equal area and equal probability mass. The partition is constructed recursively, using only lines or “fans”; the latter occurs in the lower portion of Figure 2.1c. 19 2.2 Defining x,y,α ,β ,L, andR. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23 2.3 Shifting a partition horizontally (2.3a) or vertically (2.3b). . . . . . . . . . . . . . . . . . . . 24 2.4 Service region and partitions for the uniform demand density function in the simulation study (m=16). The star represents the depot location. . . . . . . . . . . . . . . . . . . . . 28 2.5 On-time performance evaluation of region partitioning policies for the uniform demand density function (m=16,v =0.08, andq =10). . . . . . . . . . . . . . . . . . . . . . . . 31 2.6 The tested nonuniform demand density function (a) and region partitions (b)-(f). The star represents the depot location. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32 2.7 On-time performance evaluation of region partitioning policies for the nonuniform demand density function (m=12,v =0.08, andq =10). . . . . . . . . . . . . . . . . . . 33 2.8 Service region and the demand density function. The star represents the depot location and the points are customer locations identified from the data set (the x and y axes stand for the disguised longitude and latitude). . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34 2.9 Partitions of the service region in the case study (m=18). . . . . . . . . . . . . . . . . . . 35 3.1 Partition a region with different cuts . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 38 3.2 Γ -partitions forf(·), R(x) is shaded area. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39 3.3 Orthogonal cut of a squre region withk depots and uniformµ . . . . . . . . . . . . . . . . . 40 3.4 (n/2,n/2) -partition with orthogonal cut . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41 viii 3.5 Take an example ofk =5,i =1,j =4, the initial condition can be represented by case 1 or case 2. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 42 3.6 Shifting the apex inx andy direction. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 46 3.7 Convex partition and orthogonal partition. Stars indicate chosen hubs and small dots mark order locations. While the convex partition intersects roads, the orthogonal partition adheres to roads, resulting in a more logical layout. . . . . . . . . . . . . . . . . . . . . . . 49 3.8 Delivery distances (measured in kilometers) for each subregion in the case study. . . . . . 50 3.9 Convex partition and orthogonal partition at Chicago. . . . . . . . . . . . . . . . . . . . . . 53 3.10 Delivery distances (measured in kilometers) for each subregion in the case study. The subregion/hub IDs correspond to the indices in Figure 3.9. With a total of 500 orders, the travel distance for the orthogonal partition outperforms the convex partition in all subregions. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 54 ix Abstract Region partitioning plays an important role across various domains. In this thesis, we first address the challenge of optimizing region partitioning in scenarios involving a singular hub and multiple vehicles for on-time last-mile delivery. On-time last-mile delivery is expanding rapidly as people expect faster delivery of goods ranging from groceries to medicines. Managing on-time delivery systems is challenging because of the underlying uncertainties and combinatorial nature of the routing decision. In practice, the efficiency of such systems also hinges on the driver’s familiarity with the local neighborhood. This paper studies the optimal region partitioning policy to minimize the expected delivery time of customer orders in a stochastic and dynamic setting. We allow both the order locations and on-site service times to be random and generally distributed. This policy assigns every driver to a subregion, hence making sure drivers will only be dispatched to their own territories. We characterize the structure of the optimal partitioning policy and show its expected on-time performance converges to that of the flexible dispatching policy in heavy traffic. The optimal characterization features two insightful conditions that are critical to the on-time performance of last-mile delivery systems. We then develop partitioning algorithms with performance guarantees, leveraging ham sandwich cuts and 3-partitions from discrete geometry. This algorithmic development can be of independent interest for other logistics problems. We demonstrate the efficiency of the proposed region partitioning policy via numerical experiments using synthetic and real-world data sets. x Expanding our investigation, we delve into a more intricate scenario featuring multiple hubs within the service region, each equipped with an identical number of drivers. In contrast to conventional convex partitioning methods, we introduce orthogonal partitioning as an innovative alternative. By implement- ing this approach, we achieve improved alignment between partition boundaries and geographical features like road networks. Real-world case studies further attest to the effectiveness of our proposed orthogo- nal partitioning method. Through our research, we contribute to the advancement of equitable region partitioning techniques, presenting solutions for both single and multi-hub contexts. xi Chapter1 Introduction 1.1 Introduction Dividing a piece of land into sub-regions is a natural problem that belongs to many different domains within the world of operations research, such as air traffic control, congressional districting, vehicle rout- ing, facility location, urban planning, and supply chain management. Indeed, effective division of ge- ographic territory has been a fundamental societal problem since the times of antiquity. Scientifically speaking, one of the major difficulties that problems of this type pose is their intrinsically interdisciplinary nature; in order to determine an optimal partition of territory, one must combine tools from a variety of disciplines. There are various algorithms for territory districting that use computational geometry. This is distinct from, for example, explicitly discrete network-based models, such as graph partitioning and multiway cuts. The key difficulty in such problems is the balancing of allocation objectives, such as some measure of workload in the sub-regions or balanced consumption of a resource, with geometric shape constraints. As examples of shape conditions, we may require that all sub-regions be compact, simply connected (not hav- ing holes), connected, or merely measurable. Because of the extreme breadth of literature on the subject of districting, we will focus on approaches to districting problems with the following attributes: the input re- gion is a Euclidean domainR, usuallyR 2 or a subset thereof (as opposed to a more abstract space, such as 1 a graph). Second, the objective is to divideR into districts while optimizing an objective function or satis- fying a criterion. Third, the problems of interest have some practical application to districting, as opposed to being computationally oriented (e.g., finding a good mesh triangulation for numerical integration). This report contains two models. First, we study the convex single depot optimal region partition- ing policy for on-time last-mile delivery in a stochastic and dynamic setting. In contrast to the flexible dispatching policies, we partition the service region into multiple subregions and assign every driver to a subregion. Second, we extend convex single depot region partitioning policy into nonconvex multiple depots region partitioning policy. We partition the service region into multiple subregion using orthogo- nal cuts and each subregion contains a single depot. And both models assume a general demand arrival process with spatial heterogeneity, as captured by a probability density function defined on a Euclidean plane. The objective is to minimize the expected delivery time. The first model characterizes the performance of the dedicated region partitioning policy in heavy traffic, benchmarking with the flexible dispatching policy. Our analysis reveals that the performance gap between the two diminishes in heavy traffic. We then develop partitioning algorithms that are guaranteed to find the optimal region partition to any desirable error tolerance, with polynomial time complexity. Our results are further validated on synthetic and real-world data sets through simulations. When compared to alternative partitioning policies, the optimal partitioning policy delivers considerable improvement in the on-time performance (either measured with the average performance or the worst-case performance). The second model introduced a method to cut the region into sub-regions. We will analyze whether the new cut has better time complexity than existing convex partition algorithms. 1.2 RelatedLiterature We review three streams of related literature: region partitioning in logistics, stochastic and dynamic vehicle routing, and on-time delivery operations. 2 1.2.1 RegionPartitioninginLogistics Region partitioning, or districting, is a fundamental problem in logistics research, which belongs to the broad class of tactical design decisions for the management of logistics systems [15]. [29] propose a partitioning-based approach for day-to-day package delivery operations, with the goal of minimizing the total workload. [10] discuss how to partition a service region to balance the workloads of multiple vehicles following traveling salesman problem (TSP) solutions, given full information about the customer demand distribution function. When only partial distributional information is available (such as moments or empir- ical distributions), robust optimization approaches can be applied, see, e.g., [12, 11]. In a facility location setting, [13] characterize the optimal partition that minimizes the maximum workload of facilities with Voronoi diagrams. These papers only deal with the static problem and do not account for the dispatching decisions and uncertain order arrival times in a dynamic context. Recently, [3] study the region partition- ing problem for same-day delivery systems using Voronoi diagrams. They simplify the stochasticity of the problem so other tactical design issues can be addressed (including service area and fleet sizing decisions). While they focus on maximizing the number of fulfilled orders (within a day), we aim to minimize the expected delivery time by incorporating the stochastic system dynamics explicitly. The districting ideas have also been applied to shared mobility and urban retail problems [24, 8]. Notably, the existing parti- tioning literature mainly considers a partitioning problem with respect to one measurable function (e.g., travel distance), our paper is the first to balance two measurable functions in the logistics setting, which turns out to be critical to the on-time performance of last-mile delivery systems. The proposed partition- ing algorithms are guaranteed to yield optimal system performance, overcoming the limitations of Voronoi diagrams. 3 1.2.2 StochasticandDynamicVehicleRouting Our work is related to stochastic and dynamic vehicle routing problems, where order arrivals are uncertain and vehicles have to be dispatched dynamically in response to order arrivals. In terms of methodologies, our paper follows the queueing modeling framework of [5, 6]. They propose several dynamic policies based on order batching, which are proved to be asymptotically good (with constant performance guarantees). We extend their results to analyze the optimal region partitioning policy in a general setting that allows both order arrivals and service times to be spatially heterogeneous. The analysis is based on continuous approximation that seeks to replace combinatorial solution values with simpler formulas, which has been applied in other challenging routing problems [1, 9]. We leverage the famous TSP tour length formula from [4] to characterize the on-time performance of dispatching and routing policies in a tractable way. Improved dispatching and routing policies can be developed, but due to the complexity of the problem, optimality can be hardly guaranteed. For example, insertion heuristics, a priori policies, rollout methods, and approximate dynamic programming approaches have been applied to solve for dynamic dispatching and routing decisions in a stochastic environment (see, e.g., in [2, 16, 19, 21]). For the same-day delivery problem with random customer arrivals, [28] propose an approximate forward dynamic programming approach that combines scenario-based planning and a variable neighborhood search heuristic to derive real-time solutions. We refer readers to [27] for a comprehensive review of recent advances in this area. 1.2.3 On-TimeDeliveryOperations In addition to dispatching and routing, on-time delivery has been studied from various operational per- spectives. For example, in the meal delivery context, [25] study the capacity management and workforce scheduling decisions. [26] explores the dynamic pricing schemes along with routing decisions to improve the efficiency of delivery systems. When customer’s patience level is incorporated, [14] analyze the rela- tive benefit of batching as a function of market conditions in a stylized model. [20] model driver’s routing 4 behavior using machine learning and demonstrate its value for on-time delivery, which stresses the role of driver’s tacit knowledge in a practical routing environment. We note that our modeling framework is amenable to addressing the above issues. After region partitioning we can convert a multi-vehicle problem to separable single-vehicle problems, which may facilitate the modeling of these operational decisions. 1.3 ApplicationofGeometricDistricting There are many natural operational problems for which districting—from a geometric viewpoint—is a use- ful tool. One popular application is in air traffic control, as can be found in [7, 53]; here the goal is to balance the workloads within districts, interpreted as the maximum number of airplanes contained within a dis- trict at any time. The use of partitioning to balance workloads is a pervasive theme in robotics and vehicle routing as well, as seen in [20, 21, 23, 24, 42, 43]. One recurrent theme in many of these papers is the appli- cation of the foundational results in geometric probability theory, the Beardwood–Halton–Hammersley Theorem [8], which describes the length of the optimal traveling salesman tour of a set of points sampled from a probability distribution. Specifically, if {X i } is a sequence of independent samples drawn from an absolutely continuous probability distributionf(·) defined on a compact planar region R, then the length of the shortest tour through pointsX 1 ,...,X N satisfies lim N→∞ TSP({X 1 ,...,X N }) √ N =β ZZ R p f(x)dx, A consequence of this result is the fact that, if one partitionsR into districtsD 1 ,..,D N such that R D i p f(x)dx is equal for alli, then we are guaranteed that TSP tours of the samples within each district are withinO √ N Districting problems can also be used to solve facility location problems. The seminal paper contains many such examples, as do the more recent contributions. One of the key concepts is the fact that one can use a particular districting scheme—such as the Voronoi partition to reduce the set of candidate placements 5 of a facility to a finite set. Not surprisingly, computational geometry has had a significant impact on congressional districting (gerrymandering) and school district design. Here the complexity comes from the sheer multitude of factors that influence an optimal decision, such as the number of students per district (of each demographic), the commute distance between students and their assigned school, and the quality of the schools in each district. 6 Chapter2 ConvexSingle-DepotRegionPartitioning In this section, we study a restricted class of dispatching policies that rely on region partitioning. Under this policy, the service regionR is first partitioned into m subregions (districts), and then the m vehi- cles are matched to those subregions such that each subregion is assigned exactly one driver. As orders come in, each driver will be dedicated to serving orders from the assigned subregion. We call this policy the dedicated region partitioning (districting) policy (the readers shall differentiate our dedicated region partitioning policy with the partitioning policies proposed in [5], wherein their policies permit drivers serving orders from different subregions). This policy can be reckoned a divide-and-conquer dispatching algorithm. The dedicated region partitioning policy is in contrast to the dynamic vehicle routing policies that can flexibly dispatch a driver to different areas of the service region in different periods. While each driver has to be dispatched dynamically in his/her subregion (e.g., follows the single-vehicle policyπ ∗ ), the dedicated partitioning policy drops the need to coordinate different drivers. When cast as a queueing model, the comparison between the dedicated region partitioning policy and flexible dispatching policy is analogous to that between the separate (specialized) queueing system and pooling system (with complete resource pooling). 7 The dedicated region partitioning policy is often preferable for practical last-mile logistics systems due to the ease of operations and management [29, 18]. More importantly, because most on-time last- mile delivery service providers operate in urban areas that feature high-rise buildings and complicated road environments, drivers’ familiarity of the local neighborhood is critical to the system performance. The dedicated region partitioning policy facilitates the drivers’ development of tacit local knowledge and thus accelerates their performance improvement over time. In what follows, we will simply refer to this dedicated policy as region partitioning policy whenever no confusion may arise. 2.1 ProvablyGoodRegionPartitioning Our goal is to derive the optimal region partitioning policy that minimizes the expected delivery time of customer orders. Equivalently, we aim to find the optimal way of partitioning R intom subregions, so the resulting region partitioning policy yields the best on-time performance. First, we show that any region partitioning policy is optimal in light traffic, regardless of how the region is partitioned. Theorem1. Let{R 1 ,...,R m }beapartitionofR. Thecorrespondingregionpartitioningpolicyisoptimal in light traffic, i.e., the expected delivery time is the same as the optimal policy. The intuition behind the above result is that, in light traffic, the shortest possible expected delivery time is RR R ||x− x 0 ||f(x)dx+¯s, as there is no queueing and a driver can be sent immediately from the depot to serve the customer. This is the same as the expected delivery time of any region partitioning policy. In heavy traffic, the performance of a region partitioning policy hinges on how the service region is partitioned. Specifically, region partitioning decides the distribution of customer orders arriving at every subregion, which affects not only the arrival rates but also the service duration of delivery tasks in different 8 subregions. In essence, region partitioning controls the allocation of demand and capacity (service rate) concurrently through spatial bundling of delivery orders. The demand and capacity allocation will jointly determine the workload and on-time performance of the region partitioning policy. Our performance measure for a policy π in heavy traffic is T(π )/T(π ∗ ) as ρ + 2λ RR R ||x− x 0 ||f(x)dx/mvq → 1. If T(π )/T(π ∗ )→1 in heavy traffic, then we say the policy π is asymptotically optimal asπ ∗ is conjectured to be an asymptotically optimal policy in heavy traffic. Note that there are infinitely many ways to partition a region, so it is impossible to enumerate all the possible partitions. Nevertheless, we can characterize the optimal partition with three equitable conditions: Theorem2. Let{R 1 ,...,R m } be a partition ofR. If the partition satisfies ZZ R i f(x)dx= 1 m , ZZ R i p f(x)dx= 1 m ZZ R p f(x)dx, ZZ R i ||x− x 0 ||f(x)dx= 1 m ZZ R ||x− x 0 ||f(x)dx, for i ∈ {1,...,m}, then the resulting region partitioning policy is asymptotically optimal if every vehicle follows the optimal single-vehicle policy within the subregion in heavy traffic. Proof. Proof. Assume we have a partition{R 1 ,...,R m } satisfying the above three conditions, i.e, equi- table with respect to three measurable functions: f(·), p f(·), g x 0 (·), where g x 0 (x) = ||x− x 0 ||f(x). Then we employ the single-vehicle dispatching policy implied by Theorem ?? in each subregion (with 9 m = 1), i.e., each vehicle only visits areas within the assigned subregion dynamically (but the dispatch- ing is flexible within the subregion). It follows that the expected delivery time given a demand falls in subregionR i is (in heavy traffic) T i ∼ β 2 λ i (1− 1 q ) 2 h RR R i p f i (x)dx i 2 2v 2 1− ρ i − 2λ i RR R i ||x− x 0 ||f i (x)dx vq 2 , whereλ i is the arrival rate to subregioni,f i (x) is the conditional demand density function in subregion i (f i (x)=f(x)/ RR R i f(y)dy), andρ i = ¯sλ i . Note that ρ i + 2λ i RR R i ||x− x 0 ||f i (x)dx vq = ¯s λ ZZ R i f(x)dx + 2λ RR R i f(x)dx vq ZZ R i ||x− x 0 || f(x) RR R i f(y)dy dx = ¯sλ m + 2λ vq ZZ R i ||x− x 0 ||f(x)dx = ¯sλ m + 2λ mvq ZZ R ||x− x 0 ||f(x)dx=ρ + 2λ RR R ||x− x 0 ||f(x)dx mvq , where the first equality uses the fact that λ i = λ RR R i f(x)dx and the definition of f i (x), and the next two equalities follow from the first and third equitable conditions in the assumption. So it can be verified that each subregion is also in a heavy-traffic regime when the overall workload factor ρ +2λ RR R ||x− x 0 ||f(x)dx/mvq→1. 10 It follows that T i ∼ β 2 λ RR R i f(x)dx (1− 1 q ) 2 RR R i r f(x) RR R i f(y)dy dx 2 2v 2 1− ρ − 2λ RR R ||x− x 0 ||f(x)dx mvq 2 = β 2 λ m (1− 1 q ) 2 m h RR R i p f(x)dx i 2 2v 2 1− ρ − 2λ RR R ||x− x 0 ||f(x)dx mvq 2 = β 2 λ (1− 1 q ) 2 h RR R p f(x)dx i 2 2m 2 v 2 1− ρ − 2λ RR R ||x− x 0 ||f(x)dx mvq 2 , where the first two equalities are due to λ i =λ RR R i f(x)dx and the last equality follows from the second equitable condition. Therefore, the expected delivery time to the whole service region is T ∼ m X i=1 P(x∈R i )T i = β 2 λ (1− 1 q ) 2 h RR R p f(x)dx i 2 2m 2 v 2 1− ρ − 2λ RR R ||x− x 0 ||f(x)dx mvq 2 =T(π ∗ ), where we apply the law of total expectation. Because we apply the single-vehicle spatially unbiased policy π ∗ (with m = 1) in each subregion, customers within the same subregion receive the same expected waiting time. From the above analysis,T i is approximately equal toT(π ∗ ),i=1,...,m, in heavy traffic, so the expected delivery time (and waiting time) is approximately the same across different subregions. Therefore, the derived region partitioning policy is asymptotically spatially unbiased. □ The above partitioning balances the arrival rate, the local travel (TSP) distance, as well as the aver- age distance from the depot (radial travel distance) among different subregions. As such we make sure each subregion has the same workload and service rate, and thus yields the same on-time performance. However, it is generally not possible to find a convex partitioning that is equitable with respect to three 11 measurable functions inR 2 [22]. Therefore, we can not apply the above characterization to find the optimal partition. To obtain a workable characterization, we leverage the proof of Theorem 2 to show that a partition that is equitable to two measurable functions can achieve the same performance asπ ∗ and hence be optimal. Theorem3. Let{R 1 ,...,R m } be a partition ofR. If the partition satisfies ZZ R i p f(x)dx= 1 m ZZ R p f(x)dx, ZZ R i (¯svq+2||x− x 0 ||)f(x)dx= 1 m ZZ R (¯svq+2||x− x 0 ||)f(x)dx, for i ∈ {1,...,m}, then the resulting region partitioning policy is asymptotically optimal if every vehicle follows the optimal single-vehicle policy within the subregion in heavy traffic. 12 Proof. Proof. Assume we have a partition{R 1 ,...,R m } satisfying the above two conditions, i.e, equitable with respect to two measurable functions:f(·) andh x 0 (·), whereh x 0 (x)=(¯svq+2||x− x 0 ||)f(x). Then the expected delivery time given a demand falls in subregionR i is T i ∼ β 2 λ i (1− 1 q ) 2 h RR R i p f i (x)dx i 2 2v 2 1− ρ i − 2λ i RR R i ||x− x 0 ||f i (x)dx vq 2 = β 2 λ RR R i f(x)dx (1− 1 q ) 2 RR R i r f(x) RR R i f(y)dy dx 2 2v 2 1− ¯s λ RR R i f(x)dx − 2 λ RR R i f(x)dx RR R i ||x− x 0 || f(x) RR R i f(y)dy dx/vq 2 = β 2 λ (1− 1 q ) 2 h RR R i p f(x)dx i 2 2v 2 1− λ h RR R i (¯svq+2||x− x 0 ||)f(x)dx i vq ! 2 = β 2 λ (1− 1 q ) 2 h RR R p f(x)dx i 2 2m 2 v 2 1− λ [ RR R (¯svq+2||x− x 0 ||)f(x)dx] mvq 2 = β 2 λ (1− 1 q ) 2 h RR R p f(x)dx i 2 2m 2 v 2 1− ρ − 2λ RR R ||x− x 0 ||f(x)dx mvq 2 , where the second equality follows by definition (as in Theorem 2), the third equality leverages the fact that RR R i r f(x) RR R i f(y)dy dx = RR R i p f(x)dx/ q RR R i f(x)dx and RR R i ||x− x 0 || f(x) RR R i f(y)dy dx = RR R i ||x− x 0 ||f(x)dx/ RR R i f(x)dx, and the forth equality holds due to the two equitable conditions. It can be verified from the second equitable condition that each subregion falls within the heavy traffic regime whenρ +2λ RR R ||x− x 0 ||f(x)dx/mvq→1 because ρ i + 2λ i RR R i ||x− x 0 ||f i (x)dx vq = λ RR R (¯svq+2||x− x 0 ||)f(x)dx mvq =ρ + 2λ RR R ||x− x 0 ||f(x)dx mvq . Following a similar argument in the proof of Theorem 2, we can prove that this partitioning policy achieves the same performance asπ ∗ and is spatially unbiased. Specifically, because T i is approximately the same 13 asT(π ∗ ),i=1,...,m, the region partitioning policy satisfying the two equitable conditions is asymptot- ically optimal. Further, the expected waiting time is the approximately the same for all customers in the region. □ Different from the partitioning policy in Theorem 2, the new partitioning rule proposed in Theorem 3 requires neither the arrival rate nor the average distance from the depot to be equal among different subregions. Instead, it balances the whole load factor that depends on both||x− x 0 || andf(x). The in- tegral involved in the second equitable condition, RR R i (¯svq+2||x− x 0 ||)f(x)dx, represents the sum of the average round-trip travel time (between depot and the subregion) and service time per trip (multiplied by the speed). It turns out the expected delivery time would be equalized among subregions even when the arrival rates are different. Two special cases arise when either f(·) is uniform or ¯s=0. Corollary1. Iff(·) is the uniform distribution, and the partition{R 1 ,...,R m } satisfies Area(R)= 1 m Area(R), ZZ R i ||x− x 0 ||dx= 1 m ZZ R ||x− x 0 ||dx, for i ∈ {1,...,m}, then the resulting region partitioning policy is asymptotically optimal if every vehicle follows the optimal single-vehicle policy within the subregion in heavy traffic. Corollary2. If ¯s=0, and the partition{R 1 ,...,R m } satisfies ZZ R i p f(x)dx= 1 m ZZ R p f(x)dx, ZZ R i ||x− x 0 ||f(x)dx= 1 m ZZ R ||x− x 0 ||f(x)dx, for i ∈ {1,...,m}, then the resulting region partitioning policy is asymptotically optimal if every vehicle follows the optimal single-vehicle policy within the subregion in heavy traffic. 14 A critical feature of the above partitioning strategies is the inclusion of the average distance from a customer to the depot, i.e., RR R i ||x− x 0 ||f(x)dx, which is often ignored in the existing partitioning studies. This distance becomes significant for on-time delivery because it directly impacts the service rate of the delivery system, in which the drivers need to travel back and forth between customers and the depot. Remark1. [17]proposelowerboundsandupperboundsforthetraveldistanceofacapacitatedvehiclerouting problem(CVRP).BoththelowerboundandupperboundconsistofaTSPtourlengthandtheaveragedistance between a customer and the depot. Hence the above partitioning strategies also balance the CVRP bounds between subregions, which leads to similar VRP tour lengths in different subregions. According to [7], there exists a convex partition{R 1 ,...,R m } ofR that is equitable with respect to two measurable functionsµ j :R→R + ,j∈{1,2}. Hence we can prove that there always exists a convex partitioning policy that is as good asπ ∗ . Corollary3 (Existence). There exists a convex partition{R 1 ,...,R m } such that the resulting region par- titioning policy achieves the same expected delivery time as policyπ ∗ in heavy traffic. The preceding results imply that the performance of the optimal dedicated dispatching policy can achieve that of the flexible dispatching policy with complete resource pooling in heavy traffic. Hence, the benefit of pooling diminishes as demand rate goes up. The key driving force here is the average routing time required to serve a set of orders. In heavy traffic, the average routing time in a pooled system is proportional to RR R p f(x)dx (by Theorem??). In contrast, the average routing time would decrease in a dedicated system where each driver only covers a smaller areaR i . The decrease in the average routing time boosts the service rate of each driver, which compensates the loss of flexibility and resource utilization. Below we also characterize the performance ratio of the dedicated dispatching policy and flexible dis- patching policy more exactly under a special case where the demand density function is uniform and the on-site service time is zero (so ρ = 0 and σ 2 s = 0). For ease of exposition, let the average radial travel 15 distance be ¯r = RR R ||x− x 0 ||f(x)dx and the workload factor be ω = 2λ ¯r/mvq. Because ρ = 0, the heavy condition is equivalent toω→1. Proposition1. Iff(·)istheuniformdistributionandtheon-siteservicetimeiszero,theratiooftheexpected delivery time under the optimal region partitioning policy to the expected delivery time underπ ∗ satisfies (in heavy traffic) T T(π ∗ ) ≤ 1+ (m 3/2 − m)v p 2(1+qω 2 φ)(1− ω) βλ (1− 1 q ) p Area(R)+mv p 2(1+qω 2 φ)(1− ω) , whereφ>0 is a constant depending onR. The above result implies that the ratio tends to be greater for larger fleet sizes ( m). As ω → 1, the second term diminishes at a rate of √ 1− ω. The phenomenon that pooling is not always benefitial is also observed in the queueing literature where the collective forces of variability and utilization shape the performance of a queueing network. Specifi- cally, the relative advantage of a pooling system may chip away, when the variability of service duration is large [23]. However, in our case, the main underlying mechanism is based on the structure of the TSP travel distance. 2.1.1 HeterogeneousServiceTimes The former analysis is focused on the case where the service times is assumed to be homogeneous across the region. In practice, the service time distribution may be location dependent. For instance, delivering to high-rise office buildings may take longer service times than delivering to a small townhouse. In general, we can view service time as a spatial stochastic process: {S(x): x∈R}. 16 We then use¯s(x) andσ 2 s (x) to denote the mean and variance ofS(x), respectively. Again, we assume the mean and variance are finite everywhere. The original statements in [6] are built upon the i.i.d assumptions of service time and we can show that the results can be generalized to the scenario with heterogeneous and independent service times. In this case, the fraction of vehicle time spent in on-site service is ρ = λ RR R ¯s(x)f(x)dx/m. Theorem4. When service time is heterogeneous and independent across the region, within the class of spa- tially unbiased policies, T ≥ γ 2 9 λ (1+ 1 q ) 2 h RR R p f(x)dx i 2 m 2 v 2 1− ρ − 2λ RR R ||x− x 0 ||f(x)dx mvq 2 , asρ + 2λ RR R ||x− x 0 ||f(x)dx mvq →1 whereγ ≥ 2 3 √ π is a constant. Theorem 5. When service time is heterogeneous and independent across the region, there exists a provably good spatially unbiased policy ˜ π ∗ that achieves, T(˜ π ∗ )∼ β 2 λ (1− 1 q ) 2 h RR R p f(x)dx i 2 2m 2 v 2 1− ρ − 2λ RR R ||x− x 0 ||f(x)dx mvq 2 , asρ + 2λ RR R ||x− x 0 ||f(x)dx mvq →1. While Theorem 4 follows naturally from [6], the proof of Theorem 5 overcomes the dependence on the homogeneous assumption in the original proof of [6] by devising a new flexible partition-based dispatching policy, i.e., the constructive proof of the above result applies a flexible dispatching policy. Nevertheless, analogous to Theorem 3 and Corollary 3, we show that there exists a partition and corresponding dedicated region partitioning policy that can achieve the same on-time performance as flexible policy ˜ π ∗ in heavy traffic. 17 Theorem6. There exists a partition{R 1 ,...,R m } ofR that satisfies ZZ R i p f(x)dx= 1 m ZZ R p f(x)dx, ZZ R i (¯s(x)vq+2||x− x 0 ||)f(x)dx= 1 m ZZ R (¯s(x)vq+2||x− x 0 ||)f(x)dx, for i ∈ {1,...,m}, then the resulting region partitioning policy is asymptotically optimal if every vehicle follows the optimal single-vehicle policy within the subregion in heavy traffic. 2.2 OptimalPartitioningAlgorithms This section describes a computational geometric algorithm that partitionsR into districtsR i in a desirable way. Theorem 3 already establishes the key quantities of interest, and in this section, we define functions g 1 (x),g 2 (x) and their corresponding measuresµ 1 ,µ 2 as g 1 (x)= p f(x) µ 1 (Q)= ZZ Q g 1 (x)dx ¯ g 2 (x)=(¯svq+2∥x− x 0 ∥)f(x)dx µ 2 (Q)= ZZ Q g 2 (x)dx for the sake of notational convenience (when service times are heterogeneous, we can replace ¯s with the function¯s(x)). Theorem 3 states that any partition such thatµ 1 (R i )= 1 m µ 1 (R) andµ 2 (R i )= 1 m µ 2 (R) is asymptotically optimal. In this section, we describe a recursive algorithm for finding such a partition such that all districtsR i are convex, provided thatR is convex. 18 (a) (b) (c) (d) (e) Figure 2.1: A recursive partition whereR is a square, µ 1 is the density indicated by the shaded contour map, andµ 2 is uniform; the end result in 2.1e partitionsR into m = 13 pieces, all of which have equal area and equal probability mass. The partition is constructed recursively, using only lines or “fans”; the latter occurs in the lower portion of Figure 2.1c. The following Lipschitz constants will be useful in this section’s running-time analysis: define K =sup ℓ,i Z ℓ g i (x)dx, where ℓ denotes a line inR 2 , so that for any two parallel lines ℓ,ℓ ′ separated by a distance ϵ , the mass contained between them (which we will call a sliver) is at mostKϵ . We assume without loss of generality that the diameter ofR is1. The algorithm we propose is recursive, and consists of repeatedly subdividing the input regionR with either lines or “fans”, as shown in Figure 2.1. [7] give a constructive existence proof of such a partition by the following argument, assuming for notational convenience thatµ 1 (R)=µ 2 (R)=1. 19 Constructivepartitioningprocedure: 1. DivideR withm− 1 vertical lines that partitionR intom sliversS i , such thatµ 1 (S i )=1/m for alli, arranged from left to right. Define half-spaces H i = S i j=1 S i fori∈{1,...,m− 1}, so that µ i (H i )=i/m. 2. Compute signs s(i) = sign(µ 2 (H i )− µ 1 (H i )) for all i. Of course, if we have s(i) = 0 for some i∈{1,...,m− 1}, we have found an equitable partition ofR into two pieces. 3. Assumings(i) ∈ ± 1 for alli, Theorem 7 of [7] guarantees that there exists a either a pair(i ∗ ,j ∗ ) such thats(i ∗ )=s(j ∗ ) andi ∗ +j ∗ =m, or a triple(i ∗ ,j ∗ ,k ∗ ) such thats(i ∗ )=s(j ∗ )=s(k ∗ ) and i ∗ +j ∗ +k ∗ =m. 4. Theorem 5 of [7] then states that: (a) If we have a pair(i ∗ ,j ∗ ), then there exists a line that dividesR into two piecesR andL, with µ 1 (R)=µ 2 (R)=i ∗ /m µ 1 (L)=µ 2 (L)=j ∗ /m. (b) If we have a triple(i ∗ ,j ∗ ,k ∗ ), then there exists a fan, one of whose component rays is vertical and oriented downward, that dividesR into three piecesR,L, andU, such that µ 1 (R)=µ 2 (R)=i ∗ /m µ 1 (L)=µ 2 (L)=j ∗ /m µ 1 (U)=µ 2 (U)=k ∗ /m. 20 The proof above establishes our desired result that a convex equitable partition ofR intom pieces, because we can simply apply it recursively, obtaining smaller pieces until all have mass1/m of each measure. Note that if m is even, then there obviously exists a pair (i ∗ ,j ∗ ) as desired, namely (m/2,m/2); this is the famous ham sandwich theorem [beyer2004early]. Thus, we will now turn our attention to the problem of computing either a line or a fan that equitably dividesR, given a pair(i ∗ ,j ∗ ) or a triple(i ∗ ,j ∗ ,k ∗ ). 2.2.0.1 DividingR withaline This section assumes that we have a pair (i ∗ ,j ∗ ) in step 3 of the introduction of this section. Assuming thatR contains the origin (and thatR has diameter1, as mentioned earlier), we can perform a grid search forℓ ∗ over the domain(θ,r )∈[0,2π )× [0,1] via the parameterization ℓ={(x 1 ,x 2 ):x 1 cosθ +x 2 sinθ =r}. Given a pair(θ,r ), we define the division R =R∪L in the natural way, R(θ,r )={(x 1 ,x 2 ):x 1 cosθ +x 2 sinθ ≥ r} L(θ,r )={(x 1 ,x 2 ):x 1 cosθ +x 2 sinθ <r }, and define ϕ (θ,r ) as ϕ :     θ r     7→     µ 1 (R)− c 1 µ 2 (R)− c 2     where c i = i ∗ m µ i (R) is the desired mass of R; we seek a solution (θ ∗ ,r ∗ ) such that ϕ (θ ∗ ,r ∗ ) = (0,0), whose existence is guaranteed by the constructive procedure from the previous section. 21 Note that for any offset ϵ , we have ∥ϕ (θ,r )− ϕ (θ,r +ϵ )∥ 1 ≤ Kϵ by definition of K. If we have an angular offset (θ +ϵ,r ), then the difference between the two lines at an angleϵ apart is itself contained in a sliver of widthϵ (thanks to our assumption that the diameter ofR is 1), whence ∥ϕ (θ,r )− ϕ (θ +ϵ,r )∥ 1 ≤ Kϵ. LetB be a box containing the true solution(θ ∗ ,r ∗ ) with side lengthδ :=ϵ/ (2K). The preceding bounds guarantee that∥ϕ (θ,r )∥ 1 ≤ ϵ for all(θ,r )∈ B. The area of the two-dimensional search space[0,2π )× [0,1] is equal to2π . If we fill this search space with a square grid of side length δ , then we are guaranteed that some grid point lies inB. Obviously, the required number of grid pointsN isO(1/δ 2 )=O(K 2 /ϵ 2 ), so we conclude that we can compute a line that dividesµ 1 andµ 2 into the desired proportions (within a pre-specified tolerance ϵ ) usingO(K 2 /ϵ 2 ) mass computations, i.e. calls to an integration oracle. 2.2.0.2 Findingafeasibletriple This section assumes that we have a triple(i ∗ ,j ∗ ,k ∗ ) in step 3 of the introduction of this section. We seek a “fan” that dividesR into three convex piecesL,R, andU, as shown in FIGURE, where we assume that one of the rays comprising the fan is oriented downward. We can find this fan using a grid search as in the preceding section by first defining a function ϕ :R× (0,π ) 2 →[− 1,1] 4 as 22 (x, y) L R α β Figure 2.2: Defining x,y,α ,β ,L, andR. ϕ :             x y α β             7→             ϕ 1 (x,y,α,β ) ϕ 2 (x,y,α,β ) ϕ 3 (x,y,α,β ) ϕ 4 (x,y,α,β )             =             µ 1 (L)− c 1 µ 1 (R)− c 2 µ 2 (L)− d 1 µ 2 (R)− d 2             whereL andR are the left and right wedges as shown in 2.2, and c 1 = i ∗ m µ 1 (R) c 2 = k ∗ m µ 1 (R) d 1 = i ∗ m µ 2 (R) d 2 = k ∗ m µ 2 (R) are the desired masses of L and R. We seek a 4-tuple (x ∗ ,y ∗ ,α ∗ ,β ∗ ) such that f(x ∗ ,y ∗ ,α ∗ ,β ∗ ) = (0,0,0,0), i.e. an equitable partition ofR by a fan into 3 pieces (whose existence is guaranteed by the con- structive procedure on page 20). Consider the discrepancy between a given partition induced by(x,y,α,β ) and a horizontal offset (x+ϵ,y,α,β ). It is easy to see that the symetric difference between these two par- titions is contained between three slivers of width at most ϵ , as shown in 2.3a. Therefore, we see that ∥ϕ (x,y,α,β )− ϕ (x+ϵ,y,α,β )∥ 1 ≤ 3Kϵ. 23 (a) Horizontal offset (b) Vertical offset Figure 2.3: Shifting a partition horizontally (2.3a) or vertically (2.3b). It is also easy to see that a for a vertical offset (x,y+ϵ,α,β ), the difference between the two partitions is contained between two slivers of width at mostϵ , as shown in 2.3b, so that ∥ϕ (x,y,α,β )− ϕ (x,y+ϵ,α,β )∥ 1 ≤ 2Kϵ. Finally, if we have an angular offset (x,y,α + ϵ,β ), then the difference between the two partitions is obviously contained between two rays at an angleϵ apart, which is itself contained in a sliver of widthϵ thanks to our assumption of unit diameter. Therefore, ∥ϕ (x,y,α,β )− ϕ (x,y,α +ϵ,β )∥ 1 ≤ Kϵ and the same is clearly true for an offset to β as well. The search space for (x,y,α,β ) is of course equal toR× (0,π ) 2 . Because it is four-dimensional, a grid search for(x ∗ ,y ∗ ,α ∗ ,β ∗ ) will requireO(K 4 /ϵ 4 ) mass computations using brute force, although this can be improved. To see this, let (x ∗ ,y ∗ ,α ∗ ,β ∗ ) be the true optimal 4-tuple where f(x ∗ ,y ∗ ,α ∗ ,β ∗ ) = (0,0,0,0), and letB be a box about(x ∗ ,y ∗ ,α ∗ ,β ∗ ) with side length2ϵ . The preceding bounds guarantee that∥ϕ (x,y,α,β )∥ 1 ≤ 7Kϵ for all(x,y,α,β )∈B. Suppose that we discretizeR into a square grid with side length δ = ϵ/ (7K). There certainly exists at least one grid point (x i ,y i ) that lies in B for suitable values ofα andβ . Since we have fixed x i andy i , the functionsϕ 1 andϕ 3 are monotonically increasing in 24 α , and the functionsϕ 2 andϕ 4 are monotonically increasing inβ . Therefore, it is possible to narrow the search window over α and β to intervals of length δ usingO(− log 2 (δ )) = O(log 2 (K/ϵ )) iterations of bisection search forα andβ each. Hence, by performing these bisection searches over every grid point, we see that we can find a 4-tuple such that∥ϕ (x,y,α,β )∥ 1 ≤ ϵ by performingO( K 2 ϵ 2 log K ϵ ) mass calculations, i.e. calls to an integration oracle. 2.2.0.3 Recursivepartitioning The preceding two sections establish the following result: Theorem 7. Given a compact, convex regionR, probability measures µ 1 and µ 2 supported onR, and K andϵ asabove,itispossibletopartitionR intomconvexregionsR i suchthat|µ 1 (R i )− 1 m µ 1 (R)|≤ ϵ and |µ 2 (R i )− 1 m µ 2 (R)|≤ ϵ foralliusingO( K 2 m 3 ϵ 2 log Km ϵ )masscalculations,i.e. callstoanintegrationoracle ofµ 1 orµ 2 . Proof. Proof. We simply recursively partitionR using either a line (Section 2.2.0.1) or a fan (Section 2.2.0.2). We have to perform at mostm− 1 iterations of the two subroutines described in total, and each iteration requires at mostO( K 2 ϵ 2 log K ϵ ) mass calculations, each of which introduces an error of at most ϵ . If we follow this approach, then the accumulated error could be as high as mϵ . Therefore, we find it neces- sary to perform each subroutine with a desired tolerance ϵ ′ := ϵ/m , so that the entire process requires O( K 2 m 3 ϵ 2 log Km ϵ ) mass calculations as desired. □ Proposition 2. Given an ϵ -approximate equitable m-partition with respect to p f(x) and (¯svq +2||x− x 0 ||)f(x),theexpecteddeliverytimeforanorderfromR i undertheregionpartitioningpolicyT i satisfies(in heavy traffic) |T i − T(π ∗ )|≤ 4η (1− η ) 2 T(π ∗ ), ∀i=1,...,m, 25 where η = max 1 m ZZ R p f(x)dx, vq λ 1− λ RR R (¯svq+2||x− x 0 ||)f(x)dx mvq − 1 ϵ Proposition 2 implies that we can derive a dispatching policy that is sufficiently close to the optimal spatially unbiased policyπ ∗ by choosing the value ofϵ properly. 2.3 NumericalExperiments In this section we perform a set of numerical experiments to evaluate the performance of the proposed partitioning policy on both synthetic and real data sets. The partitioning algorithms and experiments were implemented in the Python programming language on an Intel i7-8700 3.20GHz Windows desktop with 16GB of RAM. We will refer to the optimal partitioning policy as the OP policy for ease of reference. 2.3.1 BenchmarkPolicies We consider the following benchmark dispatching policies for comparisons: 1. Horizontal partitioning policy. This policy applies horizontal cuts to partition the service region to m subregions. We consider two variants of it depending on the goal of partitioning: (1) Balancing demand arrival rates among subregions (HP1); (2) Balancing workload factors among subregions (HP2). Note that HP1 may not be stable when the demand traffic is high, whereas HP2 manages to be stable as long as the region-wise workload factor is less than one. Both HP1 and HP2 are not necessarily spatially unbiased. 2. Wedge partitioning policy. This policy leverages radial cuts to partition the service region, which gives rise to wedge-like subregions. Similar to HP1 and HP2, we test two different wedge partitioning 26 policies: WP1 balances demand arrival rates and WP2 balances workload factors between partitioned subregions. 3. Modified tour partitioning (MTP) policy. This is the dynamic dispatching policy ( π ∗ ) proposed in [5, 6]. MTP builds upon wedge partitions to form feasible tours of sizeq and service them FCFS using the first available vehicle. Specifically, MTP first collects batches of orders (with batch size b) from subregions defined by radial cuts so orders in each batch come from the same subregion. Then the optimal TSP tours are formed on these batches and assigned to vehicles in FCFS order. Unlike the dedicated partitioning policies, MTP allows drivers to flexibly visit different subregions in different dispatches, which serves as a theoretically lower bounding policy for comparisons. The batch sizeb is a tuning parameter for MTP. HP1, HP2, WP1, and WP2 are all dedicated partitioning policies, but they follow different partitioning rules. Because these partitioning policies are restricted to a specific type of partition (either using horizon- tal cuts or radial cuts), they can only balance one measurable function. While balancing demand arrival rates is intuitive, it can lead to unstable performance for certain subregions in heavy traffic. Therefore, we also test the case where the workload factor is balanced in partitioning (HP2 and WP2). Furthermore, we assume each driver follows the MTP policy in the assigned subregion for the studied dedicated partitioning policies. 2.3.2 SimulationStudy We conduct a simulation study to understand the performance of OP and benchmark policies under dif- ferent workload conditions. The service region is assumed to be of diamond shape with both diagonals equal to 1 (as shown in Figure 2.4a), thus having an area of0.5. The depot is located at the center(0,0) of the service region. The travel distance between two points on the service region is measured in thel 1 distance. Consequently, the points on the periphery of the service region are equidistant from the depot. 27 −0.5 0.0 0.5 −0.5 0.0 0.5 (a) Service region −0.5 0.0 0.5 −0.5 0.0 0.5 (b) OP −0.5 0.0 0.5 −0.5 0.0 0.5 (c) HP1 −0.5 0.0 0.5 −0.5 0.0 0.5 (d) HP2 Figure 2.4: Service region and partitions for the uniform demand density function in the simulation study (m=16). The star represents the depot location. The demand arrivals follow a Poisson process with rateλ . When implementing the MTP policy, we solve the optimal TSP tour using the Google OR-Tools for a batch of orders. Recall that the workload factor is defined as ρ +2λ RR R ||x− x 0 ||f(x)dx/mvq to reflect the actual traffic condition, and we vary the value ofλ to mimic different traffic conditions. 2.3.3 Uniformdensityfunction. We first set ¯s=0 andf(x) to be the uniform distribution over the service region. In this case, the optimal partition would be characterized by Corollary 1, which reduces to a partition that is equitable with respect to the area and the distance from the depot. Notably, because the service region is symmetric around the depot, the wedge partition that derives equal-area subregions using radial cuts will also equalize the distance from the depot to different subregions. This implies that the optimal partition would coincide with the wedge partitions (from WP1 and WP2), which is confirmed by our tests of the optimal partitioning algorithms (see, e.g., Figure 2.4b). The system parameters arev =0.08,q =10, andm=16 (we also test the case withm=20, of which similar results are obtained and presented in Appendix XX). The value ofλ varies from 0.1 (corresponding to a workload of 0.005) to 17.3 (corresponding to a workload of 0.9). When implementing OP, we select the value ofϵ such that the relative partition error is 1% (i.e.,ϵ =0.01× min{ RR R µ 1 (x)dx, RR R µ 2 (x)dx}/m). In this case, the running time of the proposed optimal partitioning algorithm is XXX seconds. The expected 28 delivery time is estimated from ten random simulation runs, each with a sufficiently long time span to ensure stability. Figure 2.5 summarizes the simulation results that measure the on-time performance from two perspec- tives: (a) the overall average delivery time across all subregions (T ) and (b) the longest subregion-wise aver- age delivery time (i.e., the worst-case average delivery time,max i T i ). We report the ratio of these two mea- sures to the average delivery time of MTP (T MTP ). The error bar represents the 90% confidence interval. The vertical dotted line depicts the critical workload threshold for HP, exceeding which would yield unsta- ble on-time performance. Specifically, this threshold can be computed as min i n RR R (¯svq+2||x− x 0 ||)f(x)dx/ m× RR R i (¯svq+2||x− x 0 ||)f(x)dx o , which is less than one if the workload is not balanced among subregions. In such a case (e.g., for HP1), the average delivery time to customers from certain subregions can be arbitrarily large when the traffic exceeds the threshold. As shown in Figure 2.5, the average delivery time of OP approaches that of MTP as the workload in- creases from 0.6 to 0.9, confirming the asymptotically good performance of OP. Moreover, as the workload goes to 0, the average delivery time of partitioning policies and MTP tend to be similar, which is consistent with our light traffic result. Compared to both HP1 and HP2, OP delivers shorter average delivery time and reduces the longest subregion-wise average delivery time significantly. Specifically, OP reduces the aver- age delivery time (and longest subregion-wise average delivery time) of HP1 by more than 50% when the workload is greater than 0.5. This improvement becomes more pronounced as workload increases. Com- paring the two horizontal partitioning policies, HP1 dominates HP2 in low-to-medium traffic regimes, but HP2 tends to outperform HP1 when the traffic goes high. This underlines the value of balancing workload factors that makes sure each subregion observes a stable on-time performance. The improvement of OP over HP2 varies between 12% and 42% for the on-time performance measures when the workload is greater than 0.7. Moreover, OP manages to maintain a more balanced delivery performance than HP1 and HP2 across the partitioned subregions: the longest subregion-wise average delivery time of OP is less than1.3 29 times the overall average delivery time across all workload conditions (while this ratio can be as large as 1.6 and 1.7 for HP1 and HP2, respectively). 2.3.4 Nonuniformdensityfunction. Next we test the performance of partitioning policies under the general distributed demand density func- tion and service time function. We set ¯s(x) = exp(x 1 )/3 andf(x)∝ 1/(|(x 1 +0.5)(x 1 − 0.5)|+0.1), where the demand density function is symmetric with respect to thex 2 -axis and increasing in the horizon- tal distance from the depot. The density function heatmap and partition results are shown in Figure 2.6, in which we observe that the partitioning result of OP is different than both wedge and horizontal partitions. OP adjusts the configuration of the partitioning in response to the nonuniform density function. The performance evaluation results for the five region partitioning policies are presented in Figure 2.7. There are two critical thresholds depicted in Figure 2.7: the one in (a) and (b) indicates the critical work- load threshold of HP1, while the one in (c) and (d) represents the critical workload threshold of WP1. As the workload approaches the corresponding critical thresholds, the average delivery time of certain sub- regions can be arbitrarily large under HP1 and WP1, as confirmed by Figures (2.7b) and (2.7d). In contrast, OP, HP2, and WP2 yield stable performance across all workload conditions. Among them, OP outperforms HP2 substantially across different workload conditions. Compared to WP2, OP delivers shorter average delivery time across the region and reduces the longest subregion-wise averge time in high-traffic regimes. Specifically, when the workload is greater than 0.7, OP outperforms WP2 by 2.5%-3.5% and 7-10% in terms ofT andmax i T i , respectively. The larger improvement ofmax i T i reflects the advantage of OP in short- ening the worst-case delivery time (due to the asymptotically spatial unbiasedness). We also observe that WP2 leads to better performance when the traffic is low or medium (workload between 0.2 and 0.6), and the biggest gap occurs when the workload is 0.6. We believe this observation is pertained to the chosen 30 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 Workload 1 2 3 4 5 T/T MTP OP HP1 HP2 (a) Average delivery time 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 Workload 1 2 3 4 5 6 7 8 max i T i /T MTP OP HP1 HP2 (b) Longest subregion-wise average delivery time Figure 2.5: On-time performance evaluation of region partitioning policies for the uniform demand density function (m=16,v =0.08, andq =10). 31 −0.5 0.0 0.5 −0.5 0.0 0.5 (a) Demand density function −0.5 0.0 0.5 −0.5 0.0 0.5 (b) OP −0.5 0.0 0.5 −0.5 0.0 0.5 (c) HP1 −0.5 0.0 0.5 −0.5 0.0 0.5 (d) HP2 −0.5 0.0 0.5 −0.5 0.0 0.5 (e) WP1 −0.5 0.0 0.5 −0.5 0.0 0.5 (f) WP2 Figure 2.6: The tested nonuniform demand density function (a) and region partitions (b)-(f). The star represents the depot location. 32 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 Workload 1 2 3 T/TMTP OP HP1 HP2 (a) OP, HP1, and HP2: T 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 Workload 1 2 3 4 5 6 7 maxiTi/TMTP OP HP1 HP2 (b) OP, HP1, and HP2: max i T i 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 Workload 1.0 1.5 2.0 2.5 T/TMTP OP WP1 WP2 (c) OP, WP1, and WP2: T 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 Workload 1.5 2.0 2.5 3.0 3.5 4.0 4.5 maxiTi/TMTP OP WP1 WP2 (d) OP, WP1, and WP2: max i T i Figure 2.7: On-time performance evaluation of region partitioning policies for the nonuniform demand density function (m=12,v =0.08, andq =10). demand density function. In the case study, we show that OP can outperform WP2 in medium-traffic regimes. 2.3.5 AReal-WorldCaseStudy In this subsection we conduct a case study based on a real-world last-mile delivery system. The service operator manages a physical depot in the urban business district of a major city in China. Most customer orders are placed via a smart phone app during weekdays. The company has designed a service region so only customers located within the service region are eligible for placing orders (as shown in Figure 33 −1.5 0.0 2.0 −2.0 0.0 1.5 (a) Service region −1.5 0.0 2.0 −2.0 0.0 1.5 (b) The fitted kernel density function Figure 2.8: Service region and the demand density function. The star represents the depot location and the points are customer locations identified from the data set (the x and y axes stand for the disguised longitude and latitude). 2.8a). The orders are perishable so the on-time performance is the main target for the company. More- over, because the service region comprises of high-rise office buildings and complicated local networks, a dedicated partitioning policy is preferred to yield stable and efficient delivery operations. The data set provided by the industry partner reveals the order placement time, order location (longi- tude and latitude, as coded byx 1 andx 2 ), service time, and other operational details. A disguised version of the data set (after removing sensitive information about customers) is made available by [20]. We divide the data set into a training set and a test set with a random80/20 split. The demand density function is estimated with a kernel method on the training order data. The fitted kernel density function is depicted in Figure 2.8b. The average service time function ¯s(x) is fitted using support vector regression to flexibly characterize the relationship between service time and demand location. The distance function is derived by fitting a weighted l 1 function to the actual travel distance data, i.e., the travel distance between a pair of points(x 1 ,x 2 ) is approximated asβ 1 |x 1 1 − x 1 2 |+β 2 |x 1 2 − x 2 2 | with an OLS regression (the resulting mean absolute error is 0.387 km). 34 We estimate the arrival rate λ based on the test set data. We set v = 20 km/hour according to the current practice and test different fleet compositions: m ∈ {14,16} and q = {10,15}, under which vehicles operate in medium-traffic regimes. The performance evaluation follows the same procedure of the simulation study. We test the performance of OP, WP1 and WP2, which represent the best region partitioning policies from the simulation study. The corresponding partitioning results are illustrated in Figure 2.9 form=14. We observe that the optimal partition behaves very differently than the wedge partition and horizontal partition, as it consists of subregions with varying orientations and shapes. Unlike the wedge partition and horizontal partition that only balance the arrival rates (i.e., demand density), the optimal partition has to trade-off the local travel distance and “line-haul" distance, along with the service time distribution. Table 2.1 summarizes the average delivery time for the three policies. Across all the scenarios, OP outperforms WP and HP substantially. Specifically, OP reduces the average delivery time and worst-case average deliv- ery time by 14.6% and 34.1%, respectively, over WP (the improvement over HP tends to be even higher). In particular, the average delivery time under WP and HP can exceed 140 minutes and 170 minutes for certain subregions. In contrast, OP manages to constrain the average delivery time under 100 minutes for most subregions. The significant gap between the on-time performance of OP and the benchmarks underscores the value of a careful region partitioning in practical delivery systems. −1.5 0.0 2.0 −2.0 0.0 1.5 (a) Optimal partition −1.5 0.0 2.0 −2.0 0.0 1.5 (b) Wedge partition −1.5 0.0 2.0 −2.0 0.0 1.5 (c) Horizontal partition Figure 2.9: Partitions of the service region in the case study (m=18). 35 Table 2.1: Average delivery time (in minutes) of different fleet configurations and region partitioning poli- cies in the case study. OP WP HP T max i T i T max i T i T max i T i m=14, q =20 65.97 106.62 69.62 149.46 83.30 172.79 m=16, q =15 49.26 63.40 61.91 99.04 83.89 142.98 m=18, q =10 41.40 58.80 50.66 94.52 51.44 135.62 36 Chapter3 OrthogonalPartition The motivation for this topic is based on the finding that most urban road networks are orthogonal (See Figure 3.1). In the previous section, we discussed region partitioning with straight lines and a single depot. However, in real life, most road networks are orthogonal. In order to find region partitions that follow the road connection, we want to find an optimal geometric orthogonal partition where cutting lines are orthogonal to each other and the workload of each subregion is equitable. In terms of workload optimiza- tion, we need to consider two levels: partition and routing. Based on the management level, especially in the presence of uncertainty, the service region is first partitioned into sub-regions. Then in each region, we solve each routing in subregion independently as a two-stage optimization model. In this section, we mainly focus on the first stage of optimization. We give an algorithm that takes as input a planar, sim- ply connected regionR, and probability density f(·). WithinR, there are n depot points defined as a collection P = {p 1 ,..,p n }, which are the start points of vehicle routes in sub-regions. We assume that each pointp i corresponds to exactly one vehicle. The vehicles must visit clients whose exact locations are unknown, but are assumed to be independent and identically distributed samples from densityf(·). Our goal is to partitionR into n disjoint subregions. Each vehicle is assigned to each subregion, so that the workloads in all subregions are asymptotically equal when a large number of samples are drawn. 37 (a) Partition with line (b) Partition with orthogonal lines Figure 3.1: Partition a region with different cuts 3.1 OptimalPartition As we discussed before, the problem turns out a equitable partition problem. We are given an atomic measure P for depots location and borel measure µ (·) to define TSP workload over service region. In general, we use a popular definition of equitable partition.(add citation) Definition1. LetR beacompactplanarregionandletµ (·)beanabsolutelycontinuousprobabilitydensity defined on R. LetP ={p 1 ,...,p n }⊂ R denote a set of n points. A orthogonal partition{R 1 ,...,R k } ofR into k sub-regions is said to be an equitablek-partition if we have RR R i µ (x)dA |P ∩R i | = RR R µ (x)dA k An orthogonal partition (Γ -partition) is formed either (See Figure 3.2) (i) by a horizontal line, or (ii) by a vertical line, or (iii) by a horizontal ray going to−∞ and a downward-oriented ray emanating from a common point. (iv) by a horizontal ray going to+∞ and a downward-oriented ray emanating from a common point. 38 Figure 3.2: Γ -partitions forf(·), R(x) is shaded area. If the region is partitioned with orthogonal cuts, each region will always satisfy the Definition 2. Definition 2. LetR be a compact planar region and letR i be sub-regions. Each sub-regionR i should be "relative"orthogonalconvextotheinputregionR,sothatforanytwopointsu,v∈R i ,theshortestL 1 -norm path betweenu andv is contained inR i . We are interested in partitioning a region with orthogonal cuts and each subregion is orthogonal con- vex. As shown in Figure 3.3, We have a square region with k depots andµ is uniformly distributed in the region. The region is partitioned into k orthogonal convex subregions and each subregion contains one depot and is equitable in measureµ . Next we prove that equitable partition always exists, and the partition region is relatively orthogonal convex toR. Theorem8. Givenanumbern,asimplepolygonS andprobabilitydensityµ (·)definedon S andacollection of pointsP ={p 1 ,...,p n }⊂ S where the vertices ofS and the points inP are all in general position, there exists an equitable Γ -partition, that a Γ -partition of S into n sub-regions S 1 ,...,S n with disjoint interiors, where each sub-regionS i contains exactly one point from P and satisfies R S i µ (x)dA= 1 n . We require that all points that are in a general position. If a point p i lying on the boundary of two subregionsR i andR j , it can be assigned to either region. In practice, the subregion boundaries can be 39 Figure 3.3: Orthogonal cut of a squre region withk depots and uniformµ . "perturbed" at each step of the algorithm, to ensurep i ∈ int(R i ) andp j ∈ int(R j ). Next, we consider the existence of(a,b) - equitable partition undera=b anda̸=b conditions in Lemma 2 and 8. Lemma 1. Given an even number k, a simple polygon S and probability density µ (·) defined on S and a collection of points P = {p 1 ,...,p k } ⊂ S where the vertices of S and the points in P are all in general position, there exists a( k 2 , k 2 )-balancedΓ -partition ofS into 2 sub-regionsS 1 ,S 2 . Proof. Proof. We define a region R(x) for allx∈R, such thatµ (R(x)) = 1/2 (as shown in Figure 3.4a). Let x = x 0 be the equation of the halving line for µ , then we have µ ({x|x ≥ x 0 }) = 1/2. We assign R(x 0 ) :={(z x , z y )|z x ≥ x 0 ,z y ∈ R}. For everyx < x 0 , we assignR(x) :={(z x , z y )|z x ≥ x,z y ≤ y}, wherey is chosen such thatµ (R(x))=1/2. For everyx<x 0 , we assignR(x):={(z x , z y )|z x ≥ x,z y ≥ y}, wherey is chosen such thatµ (R(x))=1/2. Whenx→−∞ ,R(x) is an horizontal line and lim x→+∞ R(x)=S/ lim x→−∞ R(x). Assume lim x→−∞ R(x))∩ P < n 2 , then we have lim x→+∞ R(x)∩P =|P|− lim x→−∞ R(x)∩P > n 2 . As we movex from−∞ to+∞, we include each point one by one. There existR(x), such thatR(x)∩P = n/2. We apply same method to prove the condition lim x→−∞ R(x))∩P > n 2 . □ 40 (a) Initial the partition with half line (b) Move the vertical line until find equitable partition Figure 3.4: (n/2,n/2) -partition with orthogonal cut In general, when k ̸= 2 m where m is an arbitrary positive integer, we can find (a,b)− equitable partition, wherea̸=b anda+b=k. Lemma2. Givenannumberk,asimplepolygonS andprobabilitydensityµ (·)definedon S andacollection of pointsP ={p 1 ,...,p n }⊂ S where the vertices ofS and the points inP are all in general position, there exists anq, such that(q,k− q)-balancedΓ -partition ofS into 2 sub-regionsS 1 ,S 2 . Proof. Proof. If k is even then a(k/2,k/2)-balancedΓ -partition exists by Lemma 1. Suppose thatk is odd and k ≥ 3. Let y = y i , i = 1,2,...,k− 1 be the equation of the horizontal line l i := {(x,y)|y = y i } such that the regionl + i above linel i have the propertyµ (l + i ) = i/k. If we can find a line that for some i, P ∩l + i = i k |P|, the partition by the linel i is(i,k− i)-balanced. We assume that, for somei the partition by l i is not (i,k− i)-balanced, which means the number of nodes in region l + i are greater or less than i k |P|. Suppose there is a pair(i,j) withi+j =k, without loss of generalityi<j, the regionsl + i andl + j follow propertiesµ (l + i )+µ (l + j )=1 andµ (l + j )= i k . There are two cases. Case1 Suppose that|l + i ∩P|>i and|l + j ∩P|>j. Hence bothl + i andl + j contain too many points relative to their areas. Their complement regionsl − i andl − j contain too few points. Similar to the proof of Lemma 1, we define a region R(x) for allx∈R, such thatµ (R(x))=i/k. Letx=x 0 be the verticle line 41 (a) Case 1 (b) Case 2 Figure 3.5: Take an example ofk = 5,i = 1,j = 4, the initial condition can be represented by case 1 or case 2. such thatµ ({x|x≥ x 0 }) = i/k. For everyx < x 0 , we assignR(x) :={(z x , z y )|z x ≥ x,z y ≤ y}, where y is chosen such thatµ (R(x)) = i/k. For everyx < x 0 , we assignR(x) :={(z x , z y )|z x ≥ x,z y ≥ y}, wherey is chosen such thatµ (R(x)) = i/k. Whenx→−∞ , lim x→+∞ R(x) = l − j andl − j ∩P < i. When x → −∞ , lim x→+∞ R(x) = l + i and l − i ∩P > i. As we move x from−∞ to +∞, we include each point one by one. There exist R(x), such that R(x)∩P = i. We apply same method to prove the condition |l + i ∩P|i and|l + j ∩P|<j. Hence bothl + i andl − j contain too many points relative to their areas. Their complement regionsl − i andl + j contain too few points. Consider a collections of lines{l c }, wherel c traverse froml i tol j . The functionϕ (l c ):=µ (l c )− 1 n |l c ∩P| is piecewise continuous, increasing on each of its components and decreasing at each discontinuity. Sinceϕ (l i )<0 andϕ (l j )>0, the intermediate value theorem guarantees the existence of a line ˜ l c where ϕ ( ˜ l c ) = 0. We apply same method to prove the condition|l + i ∩P|j. In general, we can use a similar idea to prove the casei > j, but we can always shift it to make sure i<j. □ 42 3.2 Algorithm It is easy to find an optimal partition to satisfy Definition 1 without any restriction of point measure. For this reason, we will require additional constraints of optimal partition. In our original problem assumption, we assume that our service regionR contains a finite borel measure µ (A) of measurable setA⊆ R and a point measureP = {p 1 ,...,p n } denote a set of depot points. Thus, an optimal partition should follow a constraint that each sub-regionR i contains a depot point. Beyond the restriction of two measures, we also imposed some constraints on the shape of sub-regions. At least, we should make sure for any two pointsu,v ∈R i , the shortestL1-norm path betweenu andv is contained inR i . This property is called relative orthogonal convex as Definition 2. Hence the algorithm should partition the compact planar regionR into subregionsR i to satisfy the following properties. Definition 3. Given a planar, simply connected polygonR, a point set P = {p 1 ,...,p n } ⊂ R, and a probability density µ (A) with measureable set A ⊆ R, an optimal partition ofR into n sub-regions, each sub-regionR i satisfy the following properties: (i) Each sub-regionR i contains exactly one elementp i ∈P. (ii) For any two pointsu,v∈R i , the shortestL 1 -norm path betweenu andv is contained inR i . (iii) All sub-regions satisfy RR R i µ (x)dA= 1 n RR R µ (x)dA. If we set µ (·) = p f(·), where f(·) is the client location, the partition satisfies the Theorem (BHH). In the following part, we will give the partitioning algorithm and time complexity to find the optimal partition. From Lemma 1 and 2, we already know the existence of optimal partition. In this section, we will describe an algorithm to find the equitable balanced or unbalanced 2-partition of a simple polygon with a 43 continuous measure and a point measure. The idea of the partitioning algorithm is to apply "divide-and- conquer" to partition a large region into smaller subregions iteratively. Theorem 9. Given a planar, simply connected polygonR, a point setP ={p 1 ,...,p n }⊂ R, and a prob- ability densityµ (A) with measurable setA⊆ R and the points inP are all in general position, there exists eitheranequitable(n/2,n/2)-partitionoranequitable(k,k− i)-partitionofR andP,andwecanfindsuch a partition inO(KlogK) running time. First, we have observed that orthogonal cutting will maintain the V-regions shape with the orthogonal cut based on the finding that the intersection of two V-regions is a V-region, and the intersection of the V-region and the plane is a V-region. Thus, we can always apply the orthogonal cut recursively. From Lemma 1 and 2, we know there always exists an orthogonal cutl ∗ that dividesR into two pieces and satisfies RR R 1 µ (x)dA |P ∩R 1 | = RR R 2 µ (x)dA |P ∩R 2 | = RR R µ (x)dA n (3.1) The orthogonal cut has one horizontal ray and one vertical ray, so we use an apexl ∗ = (x,y) to denote the location of the orthogonal cut. Assume thatR has diameter 1 in bothx direction andy direction, we can perform a grid search forl ∗ in the domain[0,1]× [0,1]. Given the apex(x,y), we can define the U andL under two conditions. Letx=x 0 be the halving line of the region that µ (x≥ x 0 ) = 1/2. Without loss of generality, we only consider the case thatx ≤ x 0 . We define U represents the upper polygon, andL represents the lower polygon. L(x,y)=        {(z x ,z y )|z x ≥ x, z y ≤ y} ifx≤ x 0 ; {(z x ,z y )|z x ≤ x, z y ≤ y} ifx>x 0 U(x,y)=        {(z x ,z y )|z x ≤ x, z y ≥ y} ifx≤ x 0 ; {(z x ,z y )|z x ≥ x, z y ≥ y} ifx>x 0 44 and define ϕ (x,y) as ϕ :     x y     →     µ (L)− c 1 |P ∩L|     wherec 1 =µ (R)×| P ∩L|/n, we want to find a solution such that µ (L)− c 1 =0 and|P ∩L|=k. For any offset ϵ , there is only an error term in continuous measureµ . ||µ (x,y)− µ (x,y+ϵ )||≤ Kϵ Therefore, we can calculate the Lipschitz bound ofµ . ||µ (x,y)− µ (x ′ ,y ′ )||≤ K(|x− x ′ |+|y− y ′ |) The area of two-dimensional search space [0,1]× [0,1] is 1. If we fill the two-dimensional space with square lattice points, the L1 distance between any two adjacent points approaches to p 1/N, so for any adjacent two points, the error term ofµ satisfies ||µ (x,y)− µ (x ′ ,y ′ )||≤ K 1 √ N Then the lattice should satisfyϵ ≥ K 1 √ N , then we getN =O(K/ϵ 2 ). If we apply the bisection search approach, we can reduce the running time. For fixed x, ϕ 1 (x,y) is monotonically increase iny. Then we can reduce the search space ofy to a lengths, and apply bisection search to grid optimaly in− log 2 (s). If we grid the search space[0,1] toN points, then we have the side lengths≈ 1/N for largeN, then for each grid pointx i , we can find the optimal y ∗ inO(log 2 N). 45 (a) Apex(x,y) (b) Apex offset in x andy Figure 3.6: Shifting the apex inx andy direction. From the previous analysis, in order to bound the error termϵ , the number of gridN should beK/ϵ at least. Thus, for eachx inR, we can find the optimal y ∗ inO(log 2 (K/ϵ )) in bisection search. Then we can restrict our search to a lattice grid with lengths=1/N, then the total running time isO(Klog 2 K). The algorithm can be applied recursively to partitionR with either(k/2,k/2) partition or(q,k− q) partition. 3.3 NumericalAnalysis In our earlier discussion, we introduced the concept of orthogonal partitioning as a method to divide the region that adheres to L1 convexity. This approach offers several advantages, as it closely aligns with the road network and ensures that the drivers of the hub can stay within their respective sections seamlessly. In this section, where we will conduct a series of numerical experiments to assess the effectiveness of the proposed partitioning policy. To demonstrate the effectiveness of partitioning, we evaluate the following benchmark dispatching policies for comparison. Firstly, we have the convex partitioning policy, which employs straight lines to 46 divide the region intom subregions (cite). Each of these policies ensures that the region is partitioned into m subregions, each containing one hub and possessing an equal area. One of the main purposes of partitioning is to make sure that the workload is fairly spread out among different subregions. To see how well a partitioning method is doing this, we use a measure of Traveling Salesperson Problem (TSP) tour distance. This measure helps us understand how efficient the partitioning is in terms of the distance a delivery person would need to travel. Since many road networks are set up like a grid, where streets run parallel and perpendicular to each other, we use L1 distance to calculate the distance between two points given the latitudes and longitudes. This is a practical way to estimate how far a delivery person would need to travel on a grid-like road network. To further explore this, we undertake two case studies in the cities of Los Angeles and Chicago. 3.3.1 ExperimentresultatLosAngelesCounty We first conduct a case study in Los Angeles (LA) County. Our study is based on the entire county, fea- turing multiple hubs strategically selected from various United States Postal Service (USPS) locations. To ensure efficient distribution and coverage, we adopt a specific approach for hub selection. Initially, we ran- domly choose one hub to serve as the starting point. Subsequently, we pick the next hub in a manner that maximizes the distance from all previously selected hubs. This selection process ensures that the chosen hubs are spread out effectively across the LA County region. Subsequently, we compare two different outcomes resulting from our partitioning efforts in Los Ange- les: the convex partition and the orthogonal partition. Both approaches guarantee the presence of a single hub within each subregion, while also ensuring that the regional areas are evenly distributed. You can see the visual results of these partitions in Figure 3.9. As depicted in Figure 3.7a, the convex partition’s boundaries intersect the road network extensively. This suggests that drivers would often have to travel beyond these boundaries when moving between 47 different locations. In contrast, Figure 3.7b presents the results of the orthogonal partition. In this scenario, the partition boundaries align with the road network, considerably reducing the need for drivers to venture beyond their designated regions. To assess the performance of the two partitioning methods, we conduct a test. We randomly gener- ate locations forN orders across Chicago. Subsequently, these orders are distributed among the created regions based on different partitions. Every subregion, denoted as i∈I, is associated with a unique hub labeled asi, and containsn i orders. This allocation is such that the sum of alln i values equals the total number of ordersN. We then calculate the total length of the TSP tour for each region. This tour involves a route originating from the designated hub, passing through all order locations within the region, and returning to the hub. Our evaluation operates on the premise of a singular driver allocated to each hub. This driver is tasked with efficiently navigating through the assigned region to deliver the orders. By doing so, we aim to gain insight into the effectiveness of the two partitioning strategies within practical delivery scenarios. From Figure 3.8, it becomes evident that the TSP distances for hubs from both the convex and orthog- onal partitions are approximately comparable. However, a notable distinction arises when considering the maximum distance. In this context, the convex partition exhibits a longer maximum delivery distance in hub 4 in comparison to the orthogonal partition. Furthermore, the convex partition also displays a higher degree of variance. This variability implies that the performance of the convex partition can be less stable compared to the more consistent outcomes of the orthogonal partition. Following the BHH theorem, which establishes a connection between the TSP tour distance and the region’s area, an increase in the number of orders within the service region is anticipated to result in a comparatively modest discrepancy between the TSP distances of subregions derived from the convex and orthogonal partitions. This implies that the variation in TSP tour distances between these partitions becomes less significant. Table 3.1 provides a detailed depiction of the outcomes obtained for different 48 (a) Convex Partition (b) Orthogonal Partition Figure 3.7: Convex partition and orthogonal partition. Stars indicate chosen hubs and small dots mark order locations. While the convex partition intersects roads, the orthogonal partition adheres to roads, resulting in a more logical layout. 49 Figure 3.8: Delivery distances (measured in kilometers) for each subregion in the case study. numbers of orders within the entire service region, offering a comprehensive understanding of the results. Number of Orders Convex Partition Orthogonal Partition mean(km) std mean(km) std 500 236.3 39.5 230.2 22.3 1000 321.5 34.1 313.0 13.4 2000 436.0 20.3 428.8 12.9 Table 3.1: Average TSP distance and standard deviation (std) among subregions for convex and orthogonal partitions at Los Angeles. Table 3.1 presents a comparison of the average TSP distances and their corresponding standard devia- tions (std) for subregions resulting from two different partitioning methods: the Convex Partition and the Orthogonal Partition. It shows the result under the different scenarios with varying numbers of orders in the service region. The table highlights that for all cases, the mean TSP distance and the standard devia- tion of the Orthogonal Partition are better (lower) than those of the Convex Partition. This suggests that, on average, deliveries are shorter and more consistent using the Orthogonal Partition. As the number of orders increases, the standard deviation values decrease. This implies that when there are more orders, the delivery distances become more predictable and consistent. Additionally, the table indicates that even 50 as the number of orders changes, the mean values for the Orthogonal Partition and the Convex Partition remain relatively comparable. This suggests that the Orthogonal Partition maintains its advantage over the Convex Partition in terms of average delivery distance and consistency. Overall, this table provides a clear comparison of the performance of the two partitioning methods under different scenarios of order quantities, highlighting the benefits of the Orthogonal Partition in terms of both average distance and consistency. However, as Figure 3.7, we can see that there are a lot of hilly areas in LA County with less population. Because of this, the way we divide the regions using the convex partition and orthogonal partition looks quite similar. In places with grid-like streets, the advantages of the orthogonal partition are more obvious. We proceed to conduct an experiment within the city of Chicago, which is distinguished by its grid-based urban planning. 3.3.2 ExperimentresultatChicago In this section, we conduct a case study in Chicago. Our study is based on the entire city, featuring mul- tiple hubs strategically selected from various Policy Station locations. To ensure efficient distribution and coverage, we adopt the same approach for hub selection as in LA County to ensure that the chosen hubs are spread out effectively across Chicago city. Next, we randomly select n order locations in Chicago and calculate the TSP tour based on three different partitions. We compare two distinct partitioning outcomes in Chicago, convex partition and orthogonal partition. Still, both types of partitioning ensure a single hub within each subregion while maintaining equitable regional areas. Figure 3.9 displays the results of these partitions. As seen in Figure 3.9a, the convex par- tition’s boundary intersects extensively with the road network. This implies that when drivers travel to 51 different locations, they often have to venture beyond these partition boundaries. On the other hand, Fig- ure 3.9b presents the outcome of the orthogonal partition. Here, the partition’s boundary aligns with the road network, effectively minimizing the need for drivers to travel beyond their designated regions. To assess the performance of the two partitioning methods, we conduct a test. We randomly generate locations forN orders across Chicago. From Figure 3.10 we can see the travel distance for the orthogonal partition outperforms the convex partition in all subregions. This observation emphasizes the superior performance of the orthogonal partition in real-world scenarios. This advantage stems from the orthog- onal partition’s ability to align boundaries with the road network, ensuring more efficient deliveries and optimized routes. Number of Orders Convex Partition Orthogonal Partition mean(km) std mean(km) std 500 48.5 14.5 28.6 9.2 1000 52.8 12.9 35.3 9.5 5000 76.1 11.4 61.8 7.0 Table 3.2: Average TSP distance and standard deviation (std) among subregions for convex and orthogonal partitions at Chicago. Similar to the experiment conducted in LA County, we explore different scenarios where a substantial number of orders are present within the service area. This investigation aims to validate the BHH theo- rem across two partitioning methods. The outcomes are illustrated in Table 3.2. Evidently, the orthogonal partition consistently maintains its superiority over the convex partition in terms of the average TSP dis- tance and the standard deviation values among subregions. Additionally, aligning with our expectations, the standard deviation demonstrates a decreasing trend, while the disparity between the orthogonal and convex partitions also diminishes. This outcome reinforces the effectiveness of the orthogonal partitioning approach, even in scenarios with significant order quantities. 52 (a) Convex Partition (b) Orthogonal Partition Figure 3.9: Convex partition and orthogonal partition at Chicago. 53 Figure 3.10: Delivery distances (measured in kilometers) for each subregion in the case study. The subre- gion/hub IDs correspond to the indices in Figure 3.9. With a total of 500 orders, the travel distance for the orthogonal partition outperforms the convex partition in all subregions. 54 Chapter4 Conclusions Our study introduces an orthogonal partitioning approach that closely aligns with the road network, en- hancing its practicality for real-world applications. We establish the existence of this partition and present an algorithm capable of efficiently identifying such orthogonal partitions within a polynomial time. To empirically validate the effectiveness of our proposed method, we perform numerical experiments in both LA County and Chicago. The outcomes consistently reveal that the orthogonal partitioning strategy not only leads to shorter Traveling Salesperson Problem (TSP) distances but also demonstrates greater stability compared to the convex partitioning method. This substantiates the advantages of our orthogonal partitioning approach across diverse scenarios. Currently, we’re focusing on the orthogonal partition using both an atomic measure, denoted as P , and a continuous measure called µ . We can broaden our approach to include two continuous measures, just as we do in convex partition settings. In practical terms, these measures could represent things like the straight-line distance between points and the distribution of service times across the delivery area. Furthermore, we’re assuming that the service region is straightforward and uninterrupted. There are no gaps or obstacles within it. As we progress, we can also explore scenarios where the service region contains holes or obstacles that need to be factored in. 55 Bibliography [1] Sina Ansari, Mehmet Başdere, Xiaopeng Li, Yanfeng Ouyang, and Karen Smilowitz. “Advancements in continuous approximation models for logistics and transportation systems: 1996–2016”. In: Transportation Research Part B: Methodological 107 (2018), pp. 229–252. [2] Nabila Azi, Michel Gendreau, and Jean-Yves Potvin. “A dynamic vehicle routing problem with multiple delivery routes”. In: Annals of Operations Research 199.1 (2012), pp. 103–112. [3] Dipayan Banerjee, Alan L Erera, and Alejandro Toriello. Fleet Sizing and Service Region Partitioning for Same-Day Delivery Systems. working paper. 2020. [4] Jillian Beardwood, John H Halton, and John Michael Hammersley. “The shortest path through many points”. In: Mathematical Proceedings of the Cambridge Philosophical Society. Vol. 55. 4. Cambridge University Press. 1959, pp. 299–327. [5] Dimitris J Bertsimas and Garrett Van Ryzin. “Stochastic and dynamic vehicle routing in the Euclidean plane with multiple capacitated vehicles”. In: Operations Research 41.1 (1993), pp. 60–76. [6] Dimitris J Bertsimas and Garrett Van Ryzin. “Stochastic and dynamic vehicle routing with general demand and interarrival time distributions”. In: Advances in Applied Probability 25.4 (1993), pp. 947–978. [7] Sergei Bespamyatnikh, David Kirkpatrick, and Jack Snoeyink. “Generalizing ham sandwich cuts to equitable subdivisions”. In: Discrete & Computational Geometry 24.4 (2000), pp. 605–622. [8] Junyu Cao, Chen Ma, and Wei Qi. Stall Economy: The Value of Mobility and Precision Deployment of Retail on Wheels. Available at SSRN 3711896. 2020. [9] Junyu Cao, Mariana Olvera-Cravioto, and Zuo-Jun Shen. “Last-mile shared delivery: A discrete sequential packing approach”. In: Mathematics of Operations Research 45.4 (2020), pp. 1466–1497. [10] John Gunnar Carlsson. “Dividing a territory among several vehicles”. In: INFORMS Journal on Computing 24.4 (2012), pp. 565–577. [11] John Gunnar Carlsson, Mehdi Behroozi, and Kresimir Mihic. “Wasserstein distance and the distributionally robust TSP”. In: Operations Research 66.6 (2018), pp. 1603–1624. 56 [12] John Gunnar Carlsson and Erick Delage. “Robust partitioning for stochastic multivehicle routing”. In: Operations research 61.3 (2013), pp. 727–744. [13] John Gunnar Carlsson and Raghuveer Devulapalli. “Dividing a territory among several facilities”. In: INFORMS Journal on Computing 25.4 (2013), pp. 730–742. [14] Mingliu Chen and Ming Hu. Courier dispatch in on-demand delivery. Available at SSRN. 2020. [15] Carlos F Daganzo. Logistics systems analysis. Springer Science & Business Media, 2005. [16] Justin C Goodson, Jeffrey W Ohlmann, and Barrett W Thomas. “Rollout policies for dynamic solutions to the multivehicle routing problem with stochastic demand and duration limits”. In: Operations Research 61.1 (2013), pp. 138–154. [17] Mordecai Haimovich and Alexander HG Rinnooy Kan. “Bounds and heuristics for capacitated routing problems”. In: Mathematics of operations Research 10.4 (1985), pp. 527–542. [18] Yixiao Huang, Lei Zhao, Warren B Powell, Yue Tong, and Ilya O Ryzhov. “Optimal learning for urban delivery fleet allocation”. In: Transportation Science 53.3 (2019), pp. 623–641. [19] Mathias A Klapp, Alan L Erera, and Alejandro Toriello. “The dynamic dispatch waves problem for same-day delivery”. In: European Journal of Operational Research 271.2 (2018), pp. 519–534. [20] Sheng Liu, Long He, and Zuo-Jun Max Shen. “On-time last-mile delivery: Order assignment with travel-time predictors”. In: Management Science 67.7 (2021), pp. 4095–4119. [21] Sheng Liu and Zhixing Luo. Multiperiod Dispatching and Routing for On-Time Delivery in a Dynamic and Stochastic Environment. arXiv preprint arXiv:2107.13058. 2021. [22] Chi-Yuan Lo, Jiřıé Matoušek, and William Steiger. “Algorithms for ham-sandwich cuts”. In: Discrete & Computational Geometry 11.4 (1994), pp. 433–452. [23] Avishai Mandelbaum and Martin I Reiman. “On pooling in queueing networks”. In: Management Science 44.7 (1998), pp. 971–981. [24] Wei Qi, Lefei Li, Sheng Liu, and Zuo-Jun Max Shen. “Shared mobility for last-mile delivery: Design, operational prescriptions, and environmental impact”. In: Manufacturing & Service Operations Management 20.4 (2018), pp. 737–751. [25] Damian Reyes, Alan Erera, Martin Savelsbergh, Sagar Sahasrabudhe, and Ryan O’Neil. The meal delivery routing problem. Optimization Online. 2018. [26] Marlin W Ulmer. “Dynamic pricing and routing for same-day delivery”. In: Transportation Science 54.4 (2020), pp. 1016–1033. [27] Marlin W Ulmer, Justin C Goodson, Dirk C Mattfeld, and Barrett W Thomas. “On modeling stochastic dynamic vehicle routing problems”. In: EURO Journal on Transportation and Logistics 9.2 (2020), p. 100008. 57 [28] Stacy A Voccia, Ann Melissa Campbell, and Barrett W Thomas. “The same-day delivery problem for online purchases”. In: Transportation Science 53.1 (2019), pp. 167–184. [29] Hongsheng Zhong, Randolph W Hall, and Maged Dessouky. “Territory planning and vehicle dispatching with driver learning”. In: Transportation Science 41.1 (2007), pp. 74–89. 58

Abstract (if available)

Abstract Region partitioning plays an important role across various domains. In this thesis, we first address the challenge of optimizing region partitioning in scenarios involving a singular hub and multiple vehicles for on-time last-mile delivery. On-time last-mile delivery is expanding rapidly as people expect faster delivery of goods ranging from groceries to medicines. Managing on-time delivery systems is challenging because of the underlying uncertainties and combinatorial nature of the routing decision. In practice, the efficiency of such systems also hinges on the driver's familiarity with the local neighborhood. This paper studies the optimal region partitioning policy to minimize the expected delivery time of customer orders in a stochastic and dynamic setting. We allow both the order locations and on-site service times to be random and generally distributed. This policy assigns every driver to a subregion, hence making sure drivers will only be dispatched to their own territories. We characterize the structure of the optimal partitioning policy and show its expected on-time performance converges to that of the flexible dispatching policy in heavy traffic. The optimal characterization features two insightful conditions that are critical to the on-time performance of last-mile delivery systems. We then develop partitioning algorithms with performance guarantees, leveraging ham sandwich cuts and 3-partitions from discrete geometry. This algorithmic development can be of independent interest for other logistics problems. We demonstrate the efficiency of the proposed region partitioning policy via numerical experiments using synthetic and real-world data sets.

Expanding our investigation, we delve into a more intricate scenario featuring multiple hubs within the service region, each equipped with an identical number of drivers. In contrast to conventional convex partitioning methods, we introduce orthogonal partitioning as an innovative alternative. By implementing this approach, we achieve improved alignment between partition boundaries and geographical features like road networks. Real-world case studies further attest to the effectiveness of our proposed orthogonal partitioning method. Through our research, we contribute to the advancement of equitable region partitioning techniques, presenting solutions for both single and multi-hub contexts.

Linked assets

University of Southern California Dissertations and Theses

Conceptually similar

PDF

An online cost allocation model for horizontal supply chains

PDF

Routing for ridesharing

PDF

New approaches for routing courier delivery services

PDF

Models and algorithms for the freight movement problem in drayage operations

PDF

Information design in non-atomic routing games: computation, repeated setting and experiment

PDF

Models and algorithms for pricing and routing in ride-sharing

PDF

Continuous approximation formulas for cumulative routing optimization problems

PDF

Continuous approximation formulas for location and hybrid routing/location problems

PDF

Routing problems for fuel efficient vehicles

PDF

Congestion reduction via private cooperation of new mobility services

PDF

Personalized Pareto-improving pricing-and-routing schemes with preference learning for optimum freight routing

PDF

Dynamic programming-based algorithms and heuristics for routing problems

PDF

Solving the empty container problem using double-container trucks to reduce vehicle miles

PDF

Train routing and timetabling algorithms for general networks

PDF

Integration of truck scheduling and routing with parking availability

PDF

Intelligent urban freight transportation

PDF

Optimum multimodal routing under normal condition and disruptions

PDF

The robust vehicle routing problem

PDF

Train scheduling and routing under dynamic headway control

PDF

Vehicle routing and resource allocation for health care under uncertainty

Asset Metadata

Creator Yu, Han (author)

Core Title Computational geometric partitioning for vehicle routing

School Viterbi School of Engineering

Degree Doctor of Philosophy

Degree Program Industrial and Systems Engineering

Degree Conferral Date 2023-12

Publication Date 09/11/2023

Defense Date 08/23/2023

Publisher Los Angeles, California (original), University of Southern California (original), University of Southern California. Libraries (digital)

Tag computational geometry,OAI-PMH Harvest,region partitioning,vehicle routing

Format theses (aat)

Language English

Contributor Electronically uploaded by the author (provenance)

Advisor Carlsson, John (committee chair), Dessouky, Maged (committee member), Ioannou, Petros A. (committee member)

Creator Email hyu376@usc.edu,yuhan0814@gmail.com

Permanent Link (DOI) https://doi.org/10.25549/usctheses-oUC113305012

Unique identifier UC113305012

Identifier etd-YuHan-12337.pdf (filename)

Legacy Identifier etd-YuHan-12337

Document Type Dissertation

Format theses (aat)

Rights Yu, Han

Internet Media Type application/pdf

Type texts

Source 20230911-usctheses-batch-1093 (batch), University of Southern California (contributing entity), University of Southern California Dissertations and Theses (collection)

Access Conditions The author retains rights to his/her dissertation, thesis or other graduate work according to U.S. copyright law. Electronic access is being provided by the USC Libraries in agreement with the author, as the original true and official version of the work, but does not grant the reader permission to use the work if the desired use is covered by copyright. It is the author, as rights holder, who must provide use permission if such use is covered by copyright.

Repository Name University of Southern California Digital Library

Repository Location USC Digital Library, University of Southern California, University Park Campus MC 2810, 3434 South Grand Avenue, 2nd Floor, Los Angeles, California 90089-2810, USA

Repository Email cisadmin@lib.usc.edu