9512.net

# Design of Multistar Many-to-Many Distribution Networks

Design of Multistar Many-to-Many Distribution Networks
Anton J. Kleywegt, Jinpyo Lee, and Amy R. Ward School of Industrial and Systems Engineering Georgia Institute of Technology December 31, 2006
Abstract We consider the problem of designing a distribution network to facilitate the repeated movement of shipments from many origins to many destinations. A su?cient number of the origin-destination shipments require less than the capacity of a vehicle, so that consolidation of shipments is economical. We consider the case in which consolidation takes place at terminals, and we assume each shipment moves through exactly one terminal on its way from its origin to its destination. Then, a major design decision is to determine the best number of terminals. We develop a continuous approximation method to estimate transportation costs as a function of the number of terminals. We use the continuous approximation method to choose the number of terminals that minimizes the sum of terminal cost and transportation cost. Numerical results indicate that the design resulting from the continuous approximation method facilitates operations with lower cost than those resulting from a widely used integer programming based design.

1

Introduction

In most distribution systems goods are transported from various origins to various destinations. For example, many retail chains manage distribution systems in which goods are transported from a number of suppliers to a number of retail stores. Much of this transportation takes place by truck, and it often happens that the ?ow rates for a substantial fraction of the origin-destination pairs are so small that it is not economical to send goods directly from each origin to each destination in a dedicated truck. That is, it is often more economical to consolidate the shipments of various origin-destination pairs, and transport such consolidated shipments in the same truck at the same time. There are many ways in which such consolidation can be accomplished. Next we give a number of examples. Figure 1 shows the locations of origins (circles) and destinations (squares) for a small example. Suppose that freight has to be moved from each origin to each destination in the example every week. To keep the example simple, we do not specify the volume of freight that has to be moved from each origin to each destination every week, but rather suppose that each week a truck can accommodate the total freight ?owing from up to four origins to all destinations, and the total freight ?owing from all origins to up to three destinations. The question at hand is how to design the transportation operations to move all freight each week at minimum cost. There are many alternative designs, and we give only a few examples that are qualitatively di?erent from each other. One alternative is to use out-of-vehicle consolidation, that is, to consolidate freight at terminal facilities, also called consolidation terminals, transshipment terminals, crossdocks, or transfer facilities. There are many alternative ways to design transportation operations involving terminal facilities. Here we present two simple, but important, ways, distinguished based on the number of terminals that each shipment passes through on its way from its origin to its destination. First, even if there are more than one terminal in the system, each shipment may move through only one terminal on its way from its origin to its destination. This type of system is quite widely used in the USA to distribute goods from multiple suppliers to multiple distribution centers or retail stores. In such a system it is typical for each vehicle to be based at a terminal, and to execute routes that take freight from the terminal to particular destinations, and after the freight has been delivered, to visit particular origins and collect freight at these origins and take it back to the terminal, for organization according to delivery routes and later delivery. This type of distribution network is called a multistar many-tomany distribution network, and is the type of system for which we propose a design approach in this paper. We

1

A
k 1 k 2

B
k 13

k 12 k 11 k 10

C D
k 9

k 3 k 4 k 5 k 6

E
k 7

F

k 8

Figure 1: Origins (circles) and destinations (squares) for freight ?ows.

will point out some of the advantages and disadvantages of such a system soon, but ?rst we give an example. Figure 2 shows a solution using out-of-vehicle consolidation with two terminals denoted by X and Y , in which each shipment moves through only one terminal on its way from its origin to its destination. The solution in Figure 2 consists of the following routes: X, D, B, A, 1, 2, 3, X; X, C, E, F, 6, 5, 4, X; Y, C, B, A, 13, 12, 11, 10, Y ; and Y, D, E, F, 9, 8, 7, Y . For example, consider a shipment that has to go from origin 1 to destination F . One week the shipment is picked up at origin 1 by the vehicle that executes route X, C, B, A, 1, 2, 3, X, then at terminal X the shipment is o?oaded from the vehicle and loaded onto the vehicle that will execute the route X, D, F, E, 6, 5, 4, X next, and later the shipment is delivered at destination F while the vehicle executes this route. ? A r??? ? rr ?r ¨¨ r ???? r % ¨ k k ??? 1 B B12 ?? ¨¨ rrr ?? ¨ d  s ? d z k j k 13 11 d ? d d ? ? d ? dk 2 C? ? y ? ? ? Q   d ???? Y ? ? ?' ? ?  k ?' 10 D ??? ? d  T  ??? ? T ? ? ? dE ? ? k E E E9 k EF 3 ¨  ? X ¨  ? ¨  ? ¨¨ ¨  ? k k k ¨¨ 4' 5 ?? 7r  i ¨ ?r ?? )  r k % ¨¨ ? k

6

8

Figure 2: Solution in which each shipment moves through only one terminal on its way from its origin to its destination.
Note that with such a system su?cient routes have to be chosen to ensure that all freight move from its origins to its destinations, without exceeding constraints such as vehicle capacity constraints. In Figure 2, freight can move from every origin to every destination, because each destination receives truck visits from both terminals (but each origin receives truck visits from only one terminal). In practice, it seems quite common to ensure that freight can move from every origin to every destination either by serving each origin from each terminal, or by serving each destination from each terminal. However, such a solution may lead to unnecessarily many visits

2

during which very little freight is picked up or delivered, and unnecessarily long distances traveled by vehicles. In this paper we search for more economical ways to facilitate all freight ?ows. An alternative system that uses out-of-vehicle consolidation is the following. Each origin and each destination is served from one terminal only, typically the terminal closest to the point (origin or destination). To enable freight to move from every origin to every destination, vehicles also move freight between each pair of terminals. The vehicles that move freight between the terminals are often of a di?erent type (usually larger) than the vehicles that pick up and deliver freight at origins and destinations. Such a distribution network is called a complete topology many-to-many distribution network. In such a system, each shipment moves through either one or two terminals on its way from its origin to its destination. Figure 3 shows a solution using out-of-vehicle consolidation with two terminals denoted by X and Y , in which each origin and each destination is served from one terminal only. For example, consider a shipment that has to go from origin 1 to destination F . One week the shipment is picked up at origin 1 by the vehicle that executes route X, D, B, A, 1, 2, 3, X, then at terminal X the shipment is o?oaded from the vehicle and loaded onto the vehicle that goes from terminal X to terminal Y , then at terminal Y the shipment is o?oaded from the vehicle and loaded onto the vehicle that will execute the route Y, F, 9, 8, 7, Y next, and later the shipment is delivered at destination F while the vehicle executes this route. k B12 ¨¨ rrr d  ? j k ¨ 11 X k d ? \$\$ 13 \$\$ d ? ? ? dk \$\$ 2 C\$ ? ? ? ? ? s d ? ? Y E' k \$ 10 D ? \$\$ d \$\$¨r \$\$\$ \$\$¨¨ Trr ? d % ? ? j r \$\$\$ E ¨ k ET\$\$ ' E9 k \$ 3 F  ?d X  ? d  ? d  ? k k k 4' 5 ?? d 7r  i ?r d ?? d )  ? r k ? k

1

¨ % ¨¨ k

Ar ?

r

r

B

6

8

Figure 3: Solution in which each shipment moves through only one or two terminals on its way from its origin to its destination.
There are various other systems in which each shipment moves through one or more terminals on its way from its origin to its destination. Some such systems, such as complete topology systems in which an origin or a destination may be served by more than one nearby terminals, are reviewed in Section 2. Another related design is a star topology many-to-many distribution network with a central terminal through which all freight ?ows. Speci?cally, with a star topology, each shipment travels from its origin through the pickup and delivery terminal serving the origin to the central terminal, from there to the pickup and delivery terminal serving the destination, and from there to the destination. Thus, in such a system, each shipment travels through three terminals on its way from its origin to its destination. Some practical systems, such as distribution systems for small packages, can be much more complicated than the basic system types described above, involving a hierarchy of terminals and many vehicle types traveling between terminals, with di?erent shipments moving through di?erent numbers of terminals on their way from their origins to their destinations. Complete topology, star topology, and hierarchical systems may facilitate solutions with less total travel distance than a multistar system. The transportation cost with such systems may also be less if the cost per unit freight and distance is less for the vehicles that travel between the terminals than for the pickup and delivery vehicles. However, when a shipment moves through more than one terminal on its way from its origin to its destination it requires more loading, o?oading, and additional load sorting operations, and such additional

3

handling incurs not only more cost, but also increases the risk that shipments are lost or damaged. Which of the types of systems is best depends on the importance of handling related costs relative to transportation costs. Finally, an alternative that does not require any terminals is to use in-vehicle consolidation. For example, one vehicle may pick up the loads at origins 1, 2, and 3, in that order, that have to go to destinations A, B, and D, and then drive to destinations D, B, and A, in that order, and deliver these loads. From destination A, the vehicle would next return to origin 1, ready to repeat the cycle 1, 2, 3, D, B, A, 1. Other vehicles would perform similar cycles. Su?cient cycles have to be chosen to ensure that all freight move from its origins to its destinations, without exceeding constraints such as vehicle capacity constraints. For example, if in addition to cycle 1, 2, 3, D, B, A, 1, cycle 1, 2, 3, E, F, C, 1 is performed, then all the freight originating from origins 1, 2, and 3 could be moved to its various destinations. In-vehicle consolidation is a reasonable alternative if for many of the origin-destination pairs, the amount of freight is not much less than the vehicle capacity. Figure 4 shows a solution using in-vehicle consolidation with the following cycles: 1, 2, 3, D, B, A, 1; 1, 2, 3, E, F, C, 1; 4, 5, 6, D, B, A, 4; 4, 5, 6, E, F, C, 4; 7, 8, 9, 10, A, B, D, 7; 7, 8, 9, 10, F, C, E, 7; 13, 11, 12, A, B, D, 13; and 11, 12, 13, C, E, F, 11. ' ??? r r ¨ A rr r r ??? ?r ¨? rr r j r ?? k rr % ¨ k ? rB j r 1 ?? i rr ¨12r ?r d ????  d ??? rr % ¨¨ r k E 13 k E 11 ??? d ¨ ???? d rr\$ ¨ \$\$ d 0  d k ???? ???¨\$\$\$\$ ? d r ¨ W rr 2 i ?  ???? i? ¨¨ C ??? rr ? ??   ????¨¨ d ? d ??? k ?? ??? ?  D r 10r  I d d ? ? rr ? ??? ? ? ? d o r ? r? d ? ?  ? ? ?E c rr k r ?  ? F EE k r 3 9 ? ? rr 0  0 ?? 0  rr  d   ? d ?  r c ?   rr d ? ?  ? k ' E k k rE7 4 5 ?? ?    r ?? ?  rr   q ? k j k 6 8

Figure 4: Solution using in-vehicle consolidation.
In this paper we focus on multistar many-to-many distribution systems in which each shipment moves through one terminal on its way from its origin to its destination. These systems are su?ciently widely used to justify an in-depth study focused on such systems. In addition, sometimes there is interest in determining which basic type of system is the best for a particular case, and clearly it is desirable to be able to estimate the cost of the optimal design for each type of system in order to ?nd the best overall system. Our major motivation is to develop a method to design such distribution systems. To facilitate such design, we develop a continuous approximation (CA) method to estimate the cost of the operations that will be conducted with a given design, and then we use these estimates to search for a good design. Chosen designs are then evaluated with more detailed calculations of the resulting operations. As mentioned, often a vehicle makes all the deliveries at the destinations on its route before it picks up any new loads at the origins on its route. Such a practice may increase the travel distance over what would have been possible if pickups and deliveries could have been done in any sequence. However, often this is infeasible or undesirable, sometimes because of scheduling constraints (for example, loads are delivered at the destinations in the morning and picked up from origins in the afternoon), or because loads that are picked up before all deliveries have been done would obstruct access to the loads in the back of the truck that still have to be delivered. We focus on the case in which all deliveries on a route take place before any pickups take place. Some aspects of our approach can also be used if pickups and deliveries can be done in any sequence on a route. This work was motivated by our collaboration with two companies that manage their own inbound distribution systems. One is a large retailer that obtains merchandize from many suppliers in the USA as well as overseas. For the purpose of the distribution system in the USA, the ports of import (more precisely, warehouses close to the

4

ports) are regarded as the origins of the imported freight. Thus the origins are suppliers and import warehouses in the USA. The destinations are individual retail stores. The retailer owns and operates a number of terminals. Each shipment moves through one terminal on its way from its origin to its destination with transportation provided by independent contract carriers. In the current system, each origin is served by all the terminals, and each destination is served by the terminal closest to the destination. The greater the travel distances of the carriers’ vehicles, the more the retailer has to pay the carriers, and thus the retailer has an incentive to design a distribution system that minimizes the total distance traveled, also taking into account the cost of owning and operating a terminal, and handling related costs. The retailer was considering increasing the number of terminals, and asked for design guidelines. The other company is a distributor of a large variety of products, most of which are sold to retailers, and some to industrial consumers. The origins are the distributor’s suppliers, most of whom are located in the USA and Canada, and the destinations are the distributor’s distribution centers located close to various cities in the USA. The distributor owns and operates a number of terminals, as well as a ?eet of trucks. As for the retailer, each shipment moves through one terminal on its way from its origin to its destination. Unlike the retailer, in the current system, most origins are served by one terminal, and each destination is served by all the terminals. Products are further distributed from the distributor’s distribution centers to the distributor’s customers with a separate ?eet of smaller vehicles. The ?rst challenge in designing such a system is to determine how many terminals there should be. Transportation cost is decreasing in the number of terminals, but the cost of owning and operating terminals is increasing in the number of terminals. Since we focus on systems in which each shipment moves through one terminal on its way from its origin to its destination, handling costs do not vary with the number of terminals. The second challenge in designing such a system is to determine the best locations for the chosen number of terminals. We argue that one can accurately determine the best number of terminals without exactly determining the best locations of the terminals. We also argue that to determine the best number of terminals it is important to accurately estimate how operational costs, in this case transportation costs, vary as a function of the number of terminals. As the examples above already indicate, important operational variables are the numbers of terminals that serve each origin and each destination. It will be shown that these variables have an important impact on the transportation costs, and thus should be taken into account when choosing the number of terminals. The major part of the work is the development of a new CA approach for the estimation of transportation cost, and especially linehaul cost (the average distance from a terminal to the center of the points on a route), as a function of the number of terminals as well as the numbers of terminals that serve each origin and each destination, before knowing the exact terminal locations. Then we use these estimates to search for the best number of terminals as well as the numbers of terminals that serve each origin and each destination. Thereafter the terminals are located, and the resulting design is evaluated through detailed calculations of the vehicle routes to move all the freight from its origins to its destinations. The remainder of this paper is organized as follows. We review related literature in Section 2, and we formulate the problem in Section 3. Section 4 provides a qualitative discussion of the factors in?uencing total cost, and identi?es the factors that do and do not seem important in the selection of the best number of terminals. Section 5 describes our CA method for the design of multistar many-to-many distribution networks. We also provide a procedure for making operational decisions (which terminals to use for each origin-destination ?ow and how to route the vehicles from each terminal) in Section 6 that we use in Section 7 to evaluate our solutions compared with those produced by a widely used approach.

2

Related Literature

Optimization problems that incorporate both facility location and vehicle routing decisions, such as the problem in this paper, are called location-routing problems. Laporte (1988) proposed a classi?cation of location-routing problems, and gave an overview of both exact branch-and-bound algorithms and heuristics for a number of speci?c location-routing problems. However, the paper considered only the case of a single commodity — all shipments were considered as the same product, and thus a shipment picked up from an origin did not have a

5

speci?ed destination. In the problem we consider, each shipment has a speci?ed origin-destination pair. Also, the paper did not consider the possibility of deliveries and pickups on the same route. The method that we propose is applicable both for the case with deliveries and pickups on separate routes, as well as for the case with deliveries and pickups on the same routes. The type of approach that we consider, namely a two-phase approach in which location decisions are made in the ?rst phase and routing decisions are made repeatedly, possibly with varying data, in the second phase, and future routing costs are approximated when location decisions are made, was mentioned on p. 192 as a promising research area. Laporte et al. (1989) formulated and solved integer programming models of stochastic location-routing problems. The problems involve decisions regarding the location of a single depot among a set of candidate sites, the vehicle ?eet size, and the pickup routes. All these decisions are made before the pickup quantities become known. The models constrain either the probability of a route failure or the expected penalty of a route failure due to insu?cient capacity to handle the pickups on a route. Vidal and Goetschalckx (1997) surveyed strategic production-distribution problems, with special emphasis on features that are important for international supply chains. Some of the problems in the survey involve facility location and distribution, and an overview is given of some mixed integer programming formulations and tools for solving these problems. The approach that we consider can be regarded as a continuous approximation (CA) approach. In the transportation literature, the term “continuous approximation” refers to an approach in which some problem data, typically discrete data such as the locations of origins, destinations, or facilities, are approximated with distributions, typically continuous distributions such as the uniform distribution. Kantorovitch (1942), Beckmann (1952), and Beardwood et al. (1959) did some path-breaking work on CA in transportation. Newell (1973) described the use of CA methods to provide insight into the qualitative behavior of various discrete optimization problems and sometimes to ?nd good solutions for such problems. Various researchers developed these ideas further. Here we only give a brief overview of this line of research, and we compare our work with some of the more closely related work. Langevin et al. (1996) surveyed and classi?ed the literature on CA methods in freight distribution. Daganzo (1984b) proposed a method to construct good tours in rectangles and presented formulas, similar to that of Beardwood et al. (1959), to approximate the lengths of the tours. Daganzo (1984a) proposed a cluster?rst route-second method to construct vehicle routes and presented formulas to approximate the lengths of the routes. The e?ect of the proportions of a rectangular district for each cluster on the sum of the linehaul and detour travel distances for the cluster was studied. The results in Daganzo (1984b) were used to approximate the detour travel distances. We also follow the approach of separately approximating the linehaul and detour travel distances. The approximations in Daganzo (1984a) were used to study various one-to-many distribution problems, that is, problems with one origin and many destinations. (Note that a result for a one-to-many problem gives a result for a many-to-one problem that is symmetric to it and vice versa). For example, Burns et al. (1985) developed formulas to approximate the transportation and inventory costs for one-to-many direct shipping and peddling (routes with multiple stops) distribution strategies. Daganzo (1985) derived expressions for transportation and inventory costs, argued that vehicles should be fully loaded (assuming that inventory costs are not very high), and then derived an expression for the optimal frequency with which to pick up loads for di?erent classes of items. Daganzo (1988a) showed that if products going to the di?erent destinations are homogeneous, the same vehicles operate everywhere, and there are no route length restrictions, then in-vehicle consolidation is cheaper than out-of-vehicle consolidation, that is, transshipments are not economical. Daganzo (1988b) considered a setting with more general vehicle capacity constraints than Daganzo (1988a). Also, di?erent items have di?erent characteristics, so that optimal vehicle utilization (minimizing total vehicle miles traveled) may require careful selection of the mix of items to place on di?erent vehicles. In such a setting, in contrast with Daganzo (1988a), it may be optimal to consolidate some items at a terminal into more e?cient load combinations. Campbell (1990a) considered a similar setting as Daganzo (1988a), except that Campbell (1990a) considered the case in which there are two vehicle sizes, and the larger vehicles cannot be used for delivering at the destinations. Distribution networks with transshipment terminals, in which larger vehicles transport freight from the origin directly to

6

terminals, and smaller vehicles transport freight from terminals on routes to destinations, were compared with distribution networks without transshipment terminals, in which smaller vehicles transport freight from the origin on routes to destinations. Campbell (1993b) considered the same setting as Campbell (1990a), except that in Campbell (1990a), each larger vehicle can visit only one terminal on a trip (that is, larger vehicles do not make multiple stop routes), whereas in Campbell (1993b) larger vehicles may visit multiple terminals on a route. Unlike many of the other papers, and similar to ours, Daganzo and Newell (1986) considered a problem that involves both design and operational decisions. Speci?cally, the design of a hierarchical distribution network with one origin and many destinations was studied. Also unlike many of the other papers, the case in which vehicles with di?erent transportation costs per distance operate at di?erent levels in the hierarchy was considered. As a result, it may be optimal for a shipment to undergo multiple transshipments from the origin to its destination. Continuous approximations have also been used to study many-to-many distribution problems, that is, problems with many origins and many destinations such as the problem that we consider. Daganzo (1978) considered many-to-many demand responsive transportation systems, such as taxicab systems with at most one request being transported at a time and dial-a-bus systems that allow multiple requests to be transported at a time. Three routing algorithms were modeled, and approximate expressions were derived for the average waiting times and average riding times of customers. Blumenfeld et al. (1985) derived cost expressions for the following cases: (1) Direct shipping from one origin to one destination. (2) Direct shipping from many origins to one destination. (3) Direct shipping from one origin to many destinations. (4) Direct shipping from many origins to many destinations. (5) Shipping from many origins to many destinations with all loads moving through a single consolidation terminal. Vehicle routing is not considered — shipments move directly from each origin to the terminal, and directly from the terminal to each destination. (6) Shipping from many origins to many destinations with some loads moving directly from origin to destination and other loads moving through a single consolidation terminal. Vehicle routing is not considered. Our work di?ers from Blumenfeld et al. (1985) as follows: 1. We consider the problem of shipping from many origins to many destinations only, and not the other (easier) cases. 2. In our problem, the number and locations of terminals are to be determined, but in cases (5) and (6) above, it is given that there is a single terminal with a given location. 3. We consider transportation costs and terminal costs, but not production setup costs and inventory costs. Instead, as is the case with all the applications that we have worked on, it is given that the transportation schedule repeats periodically (daily or weekly), and thus the inventory costs are not a?ected by the distribution decisions. 4. We take vehicle capacity constraints into account. 5. We consider vehicle routing — we allow a vehicle to visit multiple origins and/or multiple destinations on a route. Blumenfeld et al. (1987) described a number of simple models that were used to evaluate strategies for the distribution of parts from suppliers to General Motors assembly plants. A variety of strategies involving direct shipping, shipping through a single terminal, as well as peddling, were evaluated for the shipments from a single supplier (Delco) with 3 plants to 30 GM assembly plants. The models considered the trade-o? between transportation costs that favor fewer larger shipments and inventory costs that favor more frequent smaller shipments. A decomposition approach that exploited the small number of origins was proposed to ?nd the best combination of direct shipping and shipping through the terminal. Hall (1984), Daganzo (1987), Hall (1987), Campbell (1990b), and Campbell (1993a) are closely related to each other. These papers all considered a setting with an equal number of origins and destinations independently and uniformly distributed in a square (sometimes rectangular) region. The ?ow rate is the same for all origin-destination pairs. Terminals are arranged in a square (sometimes rectangular) grid. The following four distribution strategies were considered by several of these papers (as described in more detail below, some papers considered only some of these strategies, and some also considered other strategies):

7

1. One-terminal-nearest-terminal: Each shipment moves through exactly one terminal. All shipments from an origin move through the terminal closest to the origin irrespective of the location of the destination of the shipment. 2. One-terminal-minimum-distance: Each shipment moves through exactly one terminal. A shipment from an origin moves through one of the (typically four) terminals that surround the origin in the grid, that minimizes the travel distance from the origin through the terminal to the destination. 3. Two-terminal-nearest-terminal: Each shipment moves through one or two terminals. All shipments from an origin move through the terminal closest to the origin irrespective of the location of the destination of the shipment, and all shipments to a destination move through the terminal closest to the destination irrespective of the location of the origin of the shipment. 4. Two-terminal-minimum-distance: Each shipment moves through one or two terminals. A shipment from an origin moves through one of the (typically four) terminals that surround the origin in the grid, and a shipment to a destination moves through one of the (typically four) terminals that surround the destination in the grid, to minimize the travel distance from the origin through the terminals to the destination. Hall (1984) derived expressions for the average distance from origin to destination for the one-terminal-nearestterminal and the two-terminal-nearest-terminal strategies. The expressions were compared with average travel distances computed between the 37 largest standard metropolitan statistical areas in the United States. Daganzo (1987) considered the case of one-to-many distribution routes, without transhipment terminals, and the case of many-to-many distribution routes, without transhipment terminals, in addition to the four strategies mentioned above. Hall (1987) compared the four strategies in terms of (a) average travel distance, (b) number of terminals, and (c) number of links. For a given average travel distance, the number of terminals and number of links were regarded as measures of consolidation. Campbell (1990b) considered the case in which the transportation cost per unit distance between terminals may be less than the transportation cost per unit distance between origins and terminals and between terminals and destinations. The paper compared the two-terminal-nearest-terminal strategy, the two-terminal-minimum-distance strategy, and the two-terminal-minimum-cost strategy, that takes into account the di?erence between the local transportation cost and the transportation cost between terminals. The expressions for the minimum cost were used to derive the optimal spacing between the terminals in the grid for a given number of terminals. Unlike the other papers, the resulting service areas of the terminals were not equal. Campbell (1993a) evaluated the accuracy of the approximation formula for minimum average transportation cost derived in Campbell (1990b). The approximation formula was quite accurate, taking into account the extent with which the idealized assumptions of the approximation were violated. With two terminals, the errors were around 5%. However, as the number of terminals increased, the errors tended to increase as well. Our paper di?ers from these papers as follows: 1. These papers assume that origins and destinations are independently and uniformly distributed in a square or rectangular region, whereas we allow arbitrary distributions of origins and destinations in a rectangle. 2. These papers assumed that the number of origins equals the number of destinations, or equivalently that the density of origins equals the density of destinations, whereas we allow the numbers of origins and destinations to be di?erent. 3. These papers assumed that the ?ow rate from origins to destinations is a deterministic constant for all origin-destination pairs, whereas we allow ?ow rates that are random with di?erent marginal distributions for di?erent origin-destination pairs. 4. These papers required terminals to be located on a square or rectangular grid in the region (and thus that the number of terminals be a square number or the product of the numbers of rows and columns in the grid), whereas we allow any number of terminals, but we assume that, except for a single centrally located terminal, the terminals are independently and uniformly distributed in the rectangular region. 5. We focus on the case in which each shipment moves through exactly one terminal on its way from its origin to its destination, whereas these papers consider various strategies in which each shipment moves through

8

two terminals on its way from its origin to its destination. 6. For the case in which each shipment moves through exactly one terminal, the papers that considered this case required that each shipment moves through the terminal closest to the origin (one-terminal-nearestterminal routing), or each shipment moves through one of the (typically four) terminals closest to the origin (one-terminal-minimum-distance routing). That is, for one-terminal-nearest-terminal routing, each origin is served by vehicles from the closest terminal and each destination is served by vehicles from all terminals, and for one-terminal-minimum-distance routing, each origin is served by vehicles from the closest four terminals, and each destination is served by vehicles from all terminals. It is easily seen that it can be very ine?cient to serve each destination by vehicles from all terminals. We allow each origin and each destination to be served by vehicles from a chosen subset of terminals in addition to the central terminal, where the number of terminals that serves an origin/destination depends on the total ?ow rate from/to the origin/destination. As mentioned before, our approach allows pickups on the same routes as deliveries. Some CA work has addressed similar ideas. Hall (1991) considered a distribution problem with two terminals that serve as origins, and multiple destinations. The items originating at the two terminals are di?erent, so that items may be sent from a terminal to destinations close to the other terminal. As a result, after delivering the items originating at a terminal to some destinations, it may be better for a vehicle to next proceed to the other terminal to pick up loads there. The paper develops approximate expressions for the linehaul and detour distances if a vehicle travels from one terminal to a delivery district, and thereafter travels to the other terminal. Daganzo and Hall (1993) derived expressions to approximate the cost of vehicle routes to do both pickups and deliveries, where the deliveries on a route are completed before any pickups are done. It is assumed that the pickup and delivery points are independent and uniformly distributed in a region that can be partitioned into approximate rectangles, that each route can make at most C delivery stops, and that there is no bound on the number of pickups that a route can make. As mentioned above, very few papers address CA methods for network design, that is, to determine the number and locations of terminals. It was mentioned that Daganzo and Newell (1986) considered hierarchical network design for one-to-many distribution. In addition, Hall (1993) considered the design of a many-to-many freight distribution network in a local area, such as a metropolitan area. The area has one gateway terminal located in the center through which all shipments to and from locations outside the area pass. Unlike our problem, origins and destinations in the area are independently and uniformly distributed. Similar to our problem, pickup and delivery terminals are uniformly distributed, and the number of pickup and delivery terminals is a design variable. Each origin and each destination is served on vehicle routes from the pickup and delivery terminal in which service district it is located. Pickup and delivery routes are separate, and the vehicles used for pickup and delivery routes are di?erent from the vehicles used to transport freight between terminals. As is the case in our problem, transportation and terminal costs were considered, and the headway between successive routes was given. The optimal number of pickup and delivery terminals and the optimal number of stops on a vehicle route were determined for two distribution systems, namely a star topology and a complete topology. The costs for the two systems were compared, and a number of conclusions were made, such as that the star topology is better if the interterminal vehicles are large relative to the pickup and delivery vehicles, and if handling cost is small. Both Dasci and Verter (2001) and Ouyang and Daganzo (2006) considered problems in which a number of facilities and their locations are to be selected, and the region is to be partitioned into service areas such that each facility supplies the destinations in one of the service areas. Dasci and Verter (2001) considered only outbound transportation costs from the facilities to the destinations, in addition to ?xed facility costs and facility capacity costs. Ouyang and Daganzo (2006) considered both inbound transportation costs from a single origin to the facilities as well as outbound transportation costs from the facilities to the destinations; thus the problem considered by Ouyang and Daganzo (2006) corresponded to a special case of the one-to-many distribution network design problem considered by Daganzo and Newell (1986) in which each shipment moves through one terminal from the origin to its destination. Dasci and Verter (2001) developed a CA that required a service area size, or equivalently a terminal density, to be selected at each point in the region. Ouyang and

9

Daganzo (2006) proposed an algorithm to convert the service area size as a function of the point in the region to a solution with discrete terminal locations, and they also evaluated the di?erences between the CA costs and the costs resulting from their algorithm. In all the design problems described above, the input data include the spatial density of destination demand, and not many-to-many origin-destination ?ows as in our problem. Recently, Dasci and Laporte (2005) studied a Stackelberg game in which location decisions are made by two competitors. Each competitor can decide to locate many facilities, and the location decisions are represented with location density functions, instead of with discrete variables that specify the exact location of each facility.

3

Model Formulation

In this section we give a formulation of the problem we want to solve. The purpose is threefold. First, we want to give a precise statement of one version of the problem that we want to solve. Second, we want to determine what size instances can be solved with available software. Third, we want to introduce a common approximation to our problem that we will later compare with our CA method. Distribution operations take place repeatedly over multiple time periods. We consider distribution systems in which (1) each shipment moves through one terminal on its way from its origin to its destination, and (2) in each time period, the vehicles based at each open terminal transport goods to be delivered from the terminal to the destinations of the shipments, and thereafter pick up goods at their origins and return again to the terminals where the vehicles are based. Speci?cally, in each time period, goods have to be moved from origins in a set O to destinations in a set D. There is a ?nite set ? of ?ow scenarios. Each scenario ω ∈ ? has a probability or weight p(ω). The set ? may represent the support of an input probability distribution, or may be obtained from an input probability distribution by sampling, or may represent predictable di?erences in ?ow rates in di?erent time periods such as seasonal di?erences, or a combination of the above. For each scenario ω ∈ ?, let ?ow rate qij (ω) ≥ 0 denote the quantity of goods per time period that must be moved from origin i ∈ O to destination j ∈ D in scenario ω. Each vehicle has the same capacity Qv . There is a set XE of existing terminals, and a set XP of potential terminals. For each terminal m ∈ X := XE ∪ XP , let cm denote the di?erence in cost per time period between having terminal m open and operating, and not having terminal m open and operating. For each i, j ∈ O ∪ D ∪ X , let dij denote the cost to move a vehicle from point i to point j. We assume that the vehicle movement cost does not depend on the load carried by the vehicle. In addition to vehicle movement costs, there is a cost of Cv per time period for each vehicle based at each terminal, whether the vehicle is used or not, and a cost of cv for each vehicle that is used during a time period, independent of the distance traveled by the vehicle. We have to decide which of the existing and potential terminals should be open for all scenarios. Let binary decision variable um denote whether terminal m is open, that is, um := 1 0 if terminal m ∈ X is open otherwise. (1)

Let integer decision variable nm denote the number of vehicles assigned to terminal m. v For each scenario ω ∈ ?, we have to decide how to move each shipment from its origin to its destination through the open terminals. That is, for each origin-destination pair (i, j) ∈ O × D with qij (ω) > 0, we have to determine through which open terminal the shipment will move, and we have to decide how all the movements m will be handled with vehicle routes. Let the binary decision variable zij (ω) denote whether the shipment from i to j moves through terminal m, that is, ? ? 1 if terminal m ∈ X is used in moving the shipment in scenario ω ∈ ? m (2) zij (ω) := from origin i ∈ O to destination j ∈ D ? 0 otherwise.
m The de?nition of zij (ω) implies that the assignment of an origin-destination ?ow to a terminal may vary according to the scenario. If the assignment of origin-destination ?ows to terminals must be the same for all scenarios, m then one only needs zij variables that do not depend on ω.

10

Vehicle routing decisions are formulated in more detail later in the section. At this stage, it is su?cient to specify that τ (O , D , Q , d , nv ) denotes the optimum cost of the vehicle routing problem with a particular terminal, set O ? O of origins, set D ? D of destinations, quantities to be picked up and delivered given |O |+|D | by Q ∈ R+ (where for each origin i ∈ O , Qi denotes the quantity to be picked up at i and brought to the terminal, and for each destination j ∈ D , Qj denotes the quantity to be taken from the terminal and delivered at j), vehicle movement costs between the terminal and the considered origins and destinations given 2 by d ∈ R(1+|O |+|D |) , and nv vehicles with capacity Qv each. The arguments of interest of the function τ depend on the decision variables in (2), as follows: For each m ∈ X and z m (ω) ∈ {0, 1}|O|×|D| , let ? ? ? ? m Om (z m (ω)) := i∈O : zij (ω) > 0 ? ?
j∈D

Dm (z m (ω)) Qm (z m (ω), ω) i Qm (z m (ω), ω) j Q (z (ω), ω) d (z (ω))
m m m m

:= :=

j∈D :
i∈O

m zij (ω) > 0 m qij (ω)zij (ω)

for i ∈ Om (z m (ω)) for j ∈ Dm (z m (ω))

j∈D m (z m (ω))

:= := :=

m qij (ω)zij (ω)

i∈O m (z m (ω)) [Qm (z m (ω), ω) : l ∈ Om (z m (ω)) ∪ Dm (z m (ω))] l [dij : i, j ∈ Om (z m (ω)) ∪ Dm (z m (ω)) ∪ {m}] .

Then the optimization problem of interest is min (cm um + Cv nm ) + v
m∈X ω∈?

u∈{0,1}|X | ,nv ∈N|X |

p(ω)V (u, nv , ω)

(3)

where V (u, nv , ω) := minz(ω)
m∈X

τ (Om (z m (ω)), Dm (z m (ω)), Qm (z m (ω), ω), dm (z m (ω)), nm ) v
m zij (ω) = I{qij (ω)>0} m∈X m zij (ω) m zij (ω)

(4) (5) (6) (7)

subject to

for all i ∈ O, j ∈ D for all i ∈ O, j ∈ D, m ∈ X for all i ∈ O, j ∈ D, m ∈ X

≤ um ∈ {0, 1}

gives the minimum cost distribution plan for scenario ω given the open terminals speci?ed by u, vehicle ?eet sizes speci?ed by nv , and the origin-destination ?ows speci?ed by q(ω). Note that in the de?nition of Qm (z m (ω), ω) above, the quantities that have to be picked up and delivered in a time period in scenario ω are speci?ed by the origin-destination ?ows qij (ω). In practice, in a time period it is typical to deliver goods that were picked up in the previous time period and then to pick up goods and bring them to a terminal to be delivered in the next time period. If the origin-destination ?ows vary from week-to-week, then the total amount picked up in a particular week may not equal the total amount delivered in the same week. For example, if q is positive every second time period and zero every other period, then in alternating time periods goods will be picked up (but not delivered), and delivered (but not picked up). To capture such behavior accurately, a multistage formulation is needed instead of the two-stage formulation (3)–(7). The formulation above can accommodate various de?nitions of the function τ , and thus various types of vehicle routing constraints. Below we give a de?nition of τ that allows multiple vehicles to visit an origin or a destination during a time period, and that requires deliveries to be completed before any pickups are done on a route. This de?nition was motivated by the applications described in the introduction. Such a vehicle routing

11

problem is a combination of two problems that have been studied in the literature. The one problem is called the vehicle routing problem with backhauls; see, for example, Goetschalckx and Jacobs-Blecha (1989), Anily (1996), Potvin et al. (1996), Thangiah et al. (1996), Toth and Vigo (1997, 1999), Mingozzi et al. (1999), and Osman and Wassan (2002). The other problem is called the vehicle routing problem with split deliveries (sometimes split pickups); see, for example, Dror and Trudeau (1989), Dror et al. (1994), Mullaseril et al. (1997), Archetti et al. (2006), and Lee et al. (2006). As far as we know, the vehicle routing problem with backhauls and split pickups and deliveries has not been studied in the literature. Consider a set O of origins, set D of destinations, quantities to be picked up and delivered given by Q , vehicle movement costs between the considered origins, destinations, and a particular terminal given by d , and nv vehicles with capacity Qv each. Let 0 denote the terminal. Let V := {0} ∪ O ∪ D denote the set of nodes. Because deliveries must be completed before any pickups are done on a route, the feasible arc set on a route is A The decision variables are as follows: vk xijk aik := := ≥ 1 0 1 0 if vehicle k is in use otherwise if vehicle k travels on arc (i, j) otherwise := (i, j) ∈ (V )2 \ O × D : i = j .

0, the amount of goods picked up at i if i ∈ O or the amount of goods delivered to i if i ∈ D , by vehicle k.

Then
nv nv

τ (O , D , Q , d , nv ) :=

minx,v,a
(i,j)∈A k=1

dij xijk + cv
k=1

vk

(8) xijk for all i ∈ O ∪ D , k ∈ {1, . . . , nv } (9)

subject to
{j : (j,i)∈A}

xjik ? ≤ ?

=
{j : (i,j)∈A}

?

aik

xijk ? Qv
{j : (i,j)∈A}

for all i ∈ O ∪ D , k ∈ {1, . . . , nv } (10) for all k ∈ {1, . . . , nv } for all k ∈ {1, . . . , nv } for all i ∈ O for all j ∈ D (11) (12)

aik
i∈D

≤ Qv vk ≤ Qv vk

aik
i∈O nv

aik
k=1 nv

= Qi

(13)

ajk
k=1

= Qj xijk ≤ |S| ? 1

(14)

for all k ∈ {1, . . . , nv }, S ? O or S ? D , |S| ≥ 2 (15) (16)

{(i,j)∈A : i,j∈S}

vk xijk aik

∈ {0, 1} ∈ {0, 1} ≥ 0

for all k ∈ {1, . . . , nv }

for all (i, j) ∈ A, k ∈ {1, . . . , nv } (17) for all i ∈ O ∪ D , k ∈ {1, . . . , nv }18) (

12

With τ given by (8)–(18), problem (3)–(7) can be formulated as a two-stage mixed integer linear program, or simply as a large mixed integer linear program. Very small instances of such a problem can be solved with available software. In computational tests, instances with up to 5 origins, 5 destinations, 3 candidate terminals, and a single scenario, could be solved. Most instances with 6 origins, 6 destinations, 3 candidate terminals, and a single scenario, could not be solved — the computer’s memory was insu?cient to complete the branch-and-bound procedure. Problems in applications have hundreds of origins and destinations, and thus it seems that in the foreseeable future it will be impractical to solve the mixed integer linear program. Next we describe an approach to problem (3) that seems to be widely used. It is natural to attempt to improve the tractability of problem (3) by simplifying the vehicle routing problem (8)–(18). A popular way to do this is to model freight movements as ?ows on arcs instead of vehicle routes, and to ignore ?xed vehicle costs. The resulting problem is the following two echelon multicommodity (TEMC) location problem (for a deterministic m version, see for example Ghiani et al. 2004). Decision variable um is the same as in (1). Decision variable yij (ω) denotes the quantity of goods ?owing from origin i to destination j through terminal m in scenario ω. Then the TEMC location problem is as follows: min cm um +
m∈X ω∈?

u∈{0,1}|X |

p(ω)W (u, ω)

(19)

where W (u, ω) := miny(ω)
i∈O j∈D m∈X m (dim + dmj ) yij (ω) m yij (ω) = qij (ω) m∈X m yij (ω) m yij (ω)

(20) for all i ∈ O, j ∈ D for all i ∈ O, j ∈ D, m ∈ X for all i ∈ O, j ∈ D, m ∈ X (21) (22) (23)

subject to

≤ qij (ω)um ≥ 0

Note that, for any u ∈ {0, 1}|X | , problem (20)–(23) has an easy solution. For each origin-destination pair (i, j) ∈ O × D, let m(i, j) := arg minm∈X {dim + dmj : um = 1} denote the cheapest open (given u) terminal to use from origin i to destination j. Then W (u, ω) =
i∈O j∈D

di,m(i,j) + dm(i,j),j qij (ω)

Thus p(ω)W (u, ω)
ω∈?

=
i∈O j∈D

di,m(i,j) + dm(i,j),j
ω∈?

p(ω)qij (ω)

=
i∈O j∈D

di,m(i,j) + dm(i,j),j qij ? min
y ? i∈O j∈D m∈X

=

(dim + dmj ) yij ?m yij = qij ?m ?
m∈X yij ?m yij ?m

subject to

for all i ∈ O, j ∈ D for all i ∈ O, j ∈ D, m ∈ X for all i ∈ O, j ∈ D, m ∈ X

≤ qij um ? ≥ 0

13

where qij := ?

ω∈?

p(ω)qij (ω). Therefore problem (19)–(23) reduces to the following problem: ? ? ? ? cm um + (dim + dmj ) yij ?m min u,? y ? ?
m∈X i∈O j∈D m∈X

(24) (25) (26) (27) (28)

subject to
m∈X yij ?m yij ?m

yij = qij ?m ? ≤ qij um ? ≥ 0

for all i ∈ O, j ∈ D for all i ∈ O, j ∈ D, m ∈ X for all i ∈ O, j ∈ D, m ∈ X for all m ∈ X

um ∈ {0, 1}

Formulation (24)–(28) is used in several commercial software packages for distribution network design.

4 Qualitative Discussion of Important Factors in Distribution Network Design
Recall the ?rst stage problem (3), which can be rewritten as follows: min min
m∈X

N ∈{1,2,...} {u∈{0,1}|X | :

um =N } nv ∈N|X |

min

(cm um + Cv nm ) + v
m∈X ω∈?

p(ω)V (u, nv , ω)

(29)

= where

N ∈{1,2,...}

min

f (N )

f (N ) :=

{u∈{0,1}|X | :

min
m∈X

um =N } nv ∈N|X |

min

(cm um + Cv nm ) + v
m∈X ω∈?

p(ω)V (u, nv , ω)

(30)

In words, one can ?rst choose the number of terminals, then the locations of the terminals, and then the number of vehicles at each terminal. Of course, when one chooses the number N of terminals, one should take into account how the following optimization problems depend on N . Next we describe an approach for choosing the number N of terminals, approximating how the following optimization problems depend on N . Suppose that the locations most likely to be chosen for the terminals (typically the locations with smaller values of cm ) have approximately the same ?xed costs cm ≈ c. Then the ?rst term m∈X cm um in the objective function f can be replaced with cN . Next, note that the total number m∈X nm of vehicles needed can be estimated v quite accurately with the ?ow data and the vehicle capacity only, for example, by maxω∈? i∈O j∈D qij (ω)/Qv , so that the total ?xed vehicle cost m∈X Cv nm does not depend much on the chosen number N of terminals. v Selection of the optimal number nm of vehicles at each terminal is addressed later. Next, suppose that we approxiv mate the remaining part of the objective function, min{u∈{0,1}|X | : m∈X um =N } minnv ∈N|X | ω∈? p(ω)V (u, nv , ω), ? ? with V (N ) := p(ω)V (N, ω). Then one obtains an approximating problem
ω∈? N ∈{1,2,...}

min

? ? f (N ) := cN + V (N )

(31)

? It remains to show how an approximation V (N ) can be constructed that is both accurate and easy to compute. In the remainder of this section, we make a few observations that guide the development of our CA method regarding (1) the locations of the origins and destinations and the ?ow rates between them, (2) the importance of terminal location, and (3) the selection of how many terminals serve each origin and each destination. In most applications, the distribution of origins and destinations is quite di?erent from the uniform distribution. A relevant question is whether approximation of the locations of origins and destinations with a uniform

14

Thus the total weighted distance over all origin-destination pairs can be calculated for any given set of terminal locations. The number N of terminals was varied from 1 to 10. For each number N of terminals, the following procedure was repeated 10,000 times. One terminal was located centrally in a rectangular area covering all the origin and destination locations in our dataset. N ? 1 terminals were located independently and uniformly distributed in the rectangle that contains all origins and destinations. Then the total weighted distance over all origin-destination pairs was calculated as described above. The total cost per time period was equal to the total weighted distance in miles plus 25, 000×N for the terminals. We used 2000 units for the vehicle capacity. Table 4 shows the average total cost over the 10,000 replications, and the minimum total cost over the 10,000 replications (approximating the total cost of the best terminal locations) for each number N of terminals. The number of terminals minimizing the average cost is equal to the number of terminals minimizing the minimum cost, namely 4 terminals, which suggests that the optimum number of terminals can be determined quite accurately without optimizing the locations of the terminals. One further observation regarding terminal location is that in practice there are many factors that play a role in the selection of exact locations of facilities that do not lead to simple models, such as the local transportation infrastructure, availability of individual sites, property prices, taxes, cost of local labor, climate, aesthetics of the location, and other subjective preferences. We agree with Dasci and Laporte (2005) that models that “entail a

15

Figure 5: Origins (circles) and destinations (pluses) for one of the applications.

Number of Terminal Facilities 1 2 3 4 5 6 7 8 9 10

Mean 415,201 407,240 391,429 385,652 389,874 398,968 420,964 438,693 458,251 477,503

Cost Minimum Standard Deviation N/A N/A 317,864 48,997 305,709 55,679 300,914 53,914 320,991 50,023 342,544 44,989 359,528 47,112 379,866 44,242 401,393 39,517 418,487 36,558

Table 1: The total cost of locating N terminal facilities.

16

very precise representation of the locations” are “a ?xation in pursuit of theoretical accuracy that has relatively little practical value”. We do not attempt to model factors such as those mentioned above, and therefore we do not address the exact location of terminals. Fortunately, the optimal number of terminals does not seem to be very sensitive with respect to the exact locations to be chosen for the terminals. We also want to address the e?ect of the number of terminals that serve each origin and each destination. Recall the one-terminal-nearest-terminal distribution strategy described in Section 2. With that strategy, shipments from each origin are taken to the terminal closest to the origin, and from there the shipment is transported to its destination. Thus, with such a strategy, each origin is served by exactly one terminal, and each destination is served by all the terminals. The larger the number of terminals that serves an origin or destination, the larger the number of vehicles that have to stop at the origin or destination per time period. It was shown in Beardwood et al. (1959) that the optimal tour length increases proportionally to the square root of the number of points to be visited on the tour. Therefore, the detour distance for pickups or deliveries increases approximately proportionally to the square root of the number of terminals that serves origins or destinations. At the same time, the larger the number of terminals that serves origins or destinations, the smaller the linehaul portion of the transportation distance between origins and terminals, or terminals and destinations, can be made by careful selection of the terminal to be used to ?ow goods for each origin-destination pair. Therefore, there is a trade-o? between linehaul distance and detour distance as a function of the number of terminals that serves origins and destinations, and we want to capture this trade-o? with our approach, instead of ?xing a strategy such as the one-terminal-nearest-terminal strategy. Also, it seems reasonable that the larger the total quantity of goods that should be picked up at an origin or delivered to a destination in a time period, the larger the number of terminals that serves that origin or destination should be in the time period. As an extreme example, if the total quantity of goods that should be picked up at an origin in a time period is very small, then only one terminal should serve that origin in that time period, so that only one vehicle has to make a stop at that origin in that time period. It is typical for the total quantity of goods that should be picked up at an origin or delivered to a destination to vary signi?cantly among origins and destinations in the same time period, and also among time periods for the same origin or destination. Figure 6 shows the distribution over origins and destinations of the average total quantity of goods picked up or delivered per week for one of the applications that we worked on. Therefore, we want to allow the number of terminals that serves an origin or destination to be di?erent for di?erent origins and destinations in the same time period, and also di?erent for the same origin or destination in di?erent time periods.
1 0.9 Cumulative fraction of total flow 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 Cumulative fraction of origins and destinations

Figure 6: The distribution over origins and destinations of the average total quantity of goods picked up or delivered per week.
? Motivated by the observations above, we construct an approximation V (N, ω) that takes into account detailed

17

data regarding locations of origins and destinations, and the origin-destination ?ows, that approximate terminal locations with a uniform distribution, and that chooses the number of terminals that serve an origin or destination in a time period based on the total quantity of goods that should be picked up at the origin or delivered to the destination in the time period. Note that this approach is in some sense the opposite of the approach followed in Hall (1984), Daganzo (1987), Hall (1987), Campbell (1990b), and Campbell (1993a), where it was assumed that origins and destinations are uniformly distributed, and origin-destination ?ows are the same for all origindestination pairs, and the terminals were carefully located on a rectangular grid. Also note that the purpose of the papers above was to obtain qualitative insight without data, and not to develop a method that can produce good solutions, whereas we want to develop a method that can produce good solutions.

5

Continuous Approximation Approach

In this section we describe a CA approach for designing a distribution network in which each shipment moves through one terminal on its way from its origin to its destination. First we provide an overview of the approach in Section 5.1. Thereafter we describe the steps of the method.

5.1

Overview

It was shown in (29), (4), and (8) that we can make decisions in the following sequence, and that it should be taken into account how each decision will a?ect the subsequent decisions and objective values. The following decisions are made before the scenario ω (and thus the ?ow realization q(ω)) is known: (i) Select the number of terminals. (ii) Select the location of each terminal. (iii) Select the number of vehicles at each terminal. The following decisions are made after the ?ow realization q(ω) is known: (iv) For each origin-destination pair with ?ow qij (ω) > 0, determine which terminal will be used to move the shipment from the origin to the destination. (v) For each terminal, decide how vehicles are routed from the terminal to do the pickups and deliveries of the ?ows assigned to the terminal in the previous step. Embarking on this work, we were primarily interested in the ?rst design decision, namely how to select the best number of terminals for a distribution network. However, inspection of (3), (4), and (8), reveals that the design decision (i) should not be made without taking into account design decisions (ii) and (iii), and operational decisions (iv) and (v). It was also argued that taking design decisions (ii) and (iii) and operational decisions (iv) and (v) into account by solving mixed integer programs was impractical, and hopefully unnecessary. Thus a major part of the work is aimed at the development of tractable approximations that take design decisions (ii) and (iii) and operational decisions (iv) and (v) into account, so that design decision (i) can be made relatively easily and well. We also want to test the quality of the decisions resulting from the developed approximations. To do so, we would like to take the value of design decision (i) resulting from the approximations and that resulting from formulation (24)–(28), for each solve the mixed integer programs that produce decisions (ii)–(v), and then compare the total cost resulting from the two values of design decision (i). The challenge that the mixed integer programs that produce decisions (ii)–(v) are intractable remain, and therefore we developed heuristics to produce decisions (ii)–(v). These heuristics require signi?cantly more computational e?ort than the approximations developed to make decision (i), but are tractable enough to do the comparison for given values of decision (i). Some of these heuristics may be of interest by themselves.

18

Section 5.2 describes a CA method for selecting the number of terminals. Section 5.3 describes how the terminals are located in computational experiments to test the quality of the solutions produced with the CA method. Section 5.4 provides a heuristic for determining the vehicle ?eet sizes.

5.2

Selection of Number of Terminals

The two major cost components considered in the selection of the number of terminals are terminal ?xed cost and transportation cost. Terminal ?xed cost increases proportionally with the number of terminals. Transportation cost should decrease as the number of terminals increases, because there are more terminals to choose from when routing shipments from their origins to their destinations through terminals. Transportation cost is taken to be proportional to the transportation distance. The transportation distance results from the routes that vehicles travel from each terminal to perform the pickups and deliveries, and is partitioned into two components, namely the linehaul distance and the detour distance. The linehaul distance associated with a route is the distance from the terminal to the center of the points that are visited, and back again. The detour distance is the remaining distance on the route, and is approximately the length of a tour through the points that are visited, excluding the terminal. As explained in Section 4, as the number of terminals that serves each origin or destination increases, the average linehaul distance, or the average distance from the origin of a shipment through the terminal used for the shipment to the destination, decreases, but the total number of stops on vehicle routes increases, and thus the detour distance increases. As also explained in Section 4, the larger the total ?ow from an origin or to a destination, the larger the number of terminals that should serve that origin or destination. We control the number of terminals that serves an origin or a destination as follows. Suppose there are N open terminals. Then we select N ?1 thresholds, 0 ≤ Q1 ≤ · · · ≤ QN ?1 . Consider an origin i and a scenario ω. If j∈D qij (ω) ∈ (0, Q1 ], then 1 terminal serves origin i in scenario ω. If j∈D qij (ω) ∈ (QN ?1 , ∞), then all N terminals serve origin i in scenario ω. Otherwise, if j∈D qij (ω) ∈ (Qk?1 , Qk ], then k terminals serve origin i in scenario ω. The number of terminals that serve each destination in each scenario is determined in the same way. Although the number of terminals that serves an origin or a destination is allowed to depend on the scenario, the thresholds do not depend on the scenario. Since not all terminals serve each origin and destination, care has to be taken to ensure that for each origindestination pair (i, j) with qij (ω) > 0, there is at least one terminal that serves both i and j. We do that by designating one terminal, called the center terminal, to serve all origins and destinations. An outline of the method for selecting the number of terminals and the thresholds is as follows: 1. A method is developed to approximate total linehaul distance as a function of the number of terminals and the thresholds. This method is described in Section 5.2.1. 2. Section 5.2.2 describes a method to approximate total detour distance as a function of the number of terminals and the thresholds. 3. The approximations of total linehaul distance and total detour distance as a function of the number of terminals and the thresholds are used to search for the optimal (as measured by the approximations) number of terminals and thresholds. The search method is described in Section 5.2.3.

5.2.1

Linehaul Distance Estimation

Suppose there are N open terminals, n = 0, 1, . . . , N ? 1, located in some region that does not have to include the locations of all origins and destinations. The center terminal, used by all origins and destinations, is indexed by 0. This section describes the estimation of the total linehaul distance for given values Q := (Q1 , . . . , QN ?1 ) of the thresholds. Consider a scenario ω ∈ ?, and an origin i ∈ O. Suppose that origin i is served by a set Ni of Ni terminals, including the center terminal. The selection of the Ni ? 1 terminals that serve origin i in addition to the center terminal from the set {1, 2, . . . , N ? 1} of open terminals besides the center terminal is described later. Similarly,

19

consider a destination j ∈ D, and suppose that destination j is served by a set Nj of Nj terminals, also including the center terminal. The numbers Ni and Nj depend on the thresholds Q and on the scenario ω, but the dependence is not shown in the notation. Let Nij := Ni ∩ Nj denote the set of terminals that serve both origin i and destination j. The sets Ni , Nj , and Nij depend on the thresholds Q, on the scenario ω, and on the selection of the terminals, but as before the dependence is not shown in the notation. Let λi,n,j denote the distance from origin i through terminal n to destination j. Distance λi,n,j depends on the locations of the terminals, but the dependence is not shown in the notation. Then the minimum distance from origin i to destination j through a terminal that serves both i and j is given by Λi,j :=
n∈Nij

min λi,n,j .

(32)

Distance Λi,j depends on Nij and thus on the thresholds Q, on the scenario ω, and on the selection of the terminals; and on distances λi,n,j , and thus on the locations of the terminals. For the given number N of terminals, and the given thresholds Q, the expected total linehaul distance L(N, Q) is then calculated as follows: L(N, Q) :=
ω∈?

p(ω)
i∈O j∈D

qij (ω) E[Λi,j ] Qv

(33)

where E[Λi,j ] denotes the expected value of Λi,j with respect to random parameters involved in the selection of the terminals. Next we consider the calculation of E[Λi,j ] in greater detail. Suppose that the locations of the N ?1 terminals {1, 2, . . . , N ?1} besides the center terminal are independently and identically distributed in some region. Suppose that the set Ni \ {0} of terminals that serve origin i in addition to the center terminal are selected from terminals 1, 2, . . . , N ?1 by making Ni ?1 random selections without replacement from {1, 2, . . . , N ?1}, each time selecting each remaining element with equal probability. The set Nj \ {0} of Nj ? 1 terminals that serve destination j in addition to the center terminal are selected in the same way. This selection of terminals is independent for all origins and destinations. Note that the linehaul distance Λi,j is decreasing in both Ni and Nj , as it should be. Speci?cally, let i := ( i,1 , . . . , i,N ?1 ) and j := ( j,1 , . . . , j,N ?1 ) be two independent and identically distributed random permutations of {1, 2, . . . , N ? 1}, with distribution such that each of the (N ? 1)! permutations has probability ? ? ? 1/(N ?1)!. Consider any ni < N and nj ≤ N . Let Ni := {0, i,1 , . . . , i,ni ?1 }, Ni+ := {0, i,1 , . . . , i,ni }, Nj := ?ij := Ni ∩ Nj , N + := N + ∩ Nj , Λi,j := min ? λi,n,j , and Λ+ := min ? + λi,n,j . ? ? ? ? ? ? ? {0, j,1 , . . . , i,nj ?1 }, N ij i i,j n∈Nij n∈Nij ? ? ? ?i , Nj , Nij , and Λi,j have the same distributions as Ni , Nj , Nij , and Λi,j respectively, given that Ni = ni Then, N ? ?+ ? and Nj = nj . Also, Ni+ , Nij , and Λ+ have the same distributions as Ni , Nij , and Λi,j respectively, given that i,j ? ? ? ? ? ? Ni = ni + 1 and Nj = nj . Note that, w.p.1, Ni ? N + and Nij ? N + . Thus Λi,j ≥ Λ+ w.p.1. Hence, the
i ij i,j

conditional distribution of Λi,j given Ni = ni is stochastically decreasing in ni . < Let Nij := {n ∈ Nij : λi,n,j < λi,0,j } denote the set of terminals in Nij that give a distance from i to j strictly less than the distance from i to j through the center terminal. We calculate E[Λi,j ] by conditioning on < |Nij | and |Nij |, where for a set S, |S| denotes its cardinality:
Ni ∧Nj k?1

E[Λi,j ] =
k=1∨(Ni +Nj ?N )

P[|Nij | = k]
=0

< P |Nij | =

< |Nij | = k E Λi,j |Nij | = k, |Nij | =

(34)

(Note that it always holds that |Nij | ≥ Ni + Nj ? N , irrespective of how Ni and Nj are selected.) First, note that since the terminals in Nj \ {0} are selected without replacement and with equal probabilities from the terminals 1, 2, . . . , N ? 1, independently of Ni , it follows that the random variable |Nij | follows a hypergeometric distribution. Speci?cally, given any set Ni , the conditional probability P[|Nij | = k | Ni ] for any k ∈ {1 ∨ (Ni +

20

Nj ? N ), . . . , Ni ∧ Nj } satis?es Ni ? 1 k?1 N ? Ni Nj ? k N ?1 Nj ? 1 (Ni ? 1)!(Nj ? 1)!(N ? Ni )!(N ? Nj )! (k ? 1)!(Ni ? k)!(Nj ? k)!(N ? Ni ? Nj + k)!(N ? 1)!

P[|Nij | = k | Ni ] =

=

(It is the same as the probability of drawing k ? 1 red balls from a jar containing N ? 1 balls, Ni ? 1 of which are red, in a sample of size Nj ? 1 drawn without replacement.) Since the right side above is the same for all realizations of Ni , it follows that the probability that exactly k ∈ {1 ∨ (Ni + Nj ? N ), . . . , Ni ∧ Nj } terminals serve both i and j is P[|Nij | = k] = (Ni ? 1)!(Nj ? 1)!(N ? Ni )!(N ? Nj )! (k ? 1)!(Ni ? k)!(Nj ? k)!(N ? Ni ? Nj + k)!(N ? 1)!

As before, P[|Nij | = k] depends on the thresholds Q and the scenario ω. Next, since the locations of terminals {1, 2, . . . , N ? 1} are independent and identically distributed, and the terminals in Ni \ {0} and Nj \ {0} are selected independently of the locations of the terminals, it follows that given that |Nij | = k, the number of the k ? 1 terminals n other than the center terminal that have distance λi,n,j strictly less than λi,0,j has a binomial distribution. Speci?cally, for ∈ {0, 1, . . . , k ? 1},
< P |Nij | =

|Nij | = k

=

k?1

P[λi,1,j < λi,0,j ] P[λi,1,j ≥ λi,0,j ]k?1? ∈ {0, 1, . . . , k ? 1} that

Since Λi,j is a nonnegative random variable, it holds for
< E Λi,j |Nij | = k, |Nij | = ∞

=
0 ∞

< P Λi,j > α |Nij | = k, |Nij | =

dα dα
,j

=
0

P

n∈Nij

< min λi,n,j > α |Nij | = k, |Nij | =

λi,0,j

=
0 λi,0,j

P min{λi,1,j , . . . , λi, ,j } > α λi,1,j < λi,0,j , . . . , λi, P λi,1,j > α, . . . , λi,
,j

< λi,0,j dα

=
0 λi,0,j

> α λi,1,j < λi,0,j , . . . , λi, dα

,j

< λi,0,j dα

=
0

P λi,1,j > α λi,1,j < λi,0,j

The third and ?fth equalities follow from terminals {1, 2, . . . , N ?1} being independent and identically distributed, and the terminals in Ni \ {0} and Nj \ {0} being selected independently of the locations of the terminals. As a special case, note that it holds for = 0 that
< E Λi,j |Nij | = k, |Nij | = 0

= λi,0,j

In summary, it follows from (34) that
Ni ∧Nj

E[Λi,j ] =
k=1∨(Ni +Nj ?N )

(Ni ? 1)!(Nj ? 1)!(N ? Ni )!(N ? Nj )! (Ni ? k)!(Nj ? k)!(N ? Ni ? Nj + k)!(N ? 1)!
k?1

×
=0

1 P[λi,1,j < λi,0,j ] P[λi,1,j ≥ λi,0,j ]k?1? !(k ? 1 ? )! ×
λi,0,j 0

P λi,1,j > α λi,1,j < λi,0,j

(35)

21

It remains to explain how P[λi,1,j < λi,0,j ] and P [λi,1,j > α | λi,1,j < λi,0,j ] for α ∈ (0, λi,0,j ) can be computed. We provide the details of these calculations in the appendix, where it is shown that if the terminals {1, 2, . . . , N ? 1} are uniformly distributed in a rectangular region, and distance λi,1,j is given by the L1 -metric, then P [λi,1,j > α | λi,1,j < λi,0,j ] is a piecewise polynomial in α of degree at most two, and thus
λi,0,j 0

P λi,1,j > α λi,1,j < λi,0,j

is conceptually simple. Thus, for a given number N of terminals, and given thresholds Q, the expected total linehaul distance L(N, Q) can be calculated quite easily using (33).

5.2.2

Detour Distance Estimation

Beardwood et al. (1959) established the following result that is widely used in the CA literature. Consider an independent sequence of uniformly distributed points in a set A ? Rk with Lebesque measure ?(A) > 0. Let T n denote the shortest tour length, as measured by the L2 Euclidean distance, through the ?rst n points in the sequence. Then there exists a constant βk , independent of the sequence and of A, such that with probability 1, limn→∞ n?(k?1)/k T n = βk k 1/2 [?(A)]1/k . Speci?cally, for R2 , there exists a constant β2 ∈ [0.44, 0.65] (the exact value is not yet known), such that with probability 1, limn→∞ n?1/2 T n = β2 21/2 [?(A)]1/2 . The approximation T n ≈ β n?(A), for some β depending on the distance metric, has been used in many vehicle routing applications, and in this section we do the same for the purpose of obtaining a tractable approximation of the detour distance as a function of the number N of terminals and the thresholds Q. Consider a given number N of terminals, given values Q := (Q1 , . . . , QN ?1 ) of the thresholds, and a given scenario ω ∈ ?. Assume that the vehicles are fully loaded with goods to be delivered when departing from a terminal and fully loaded with goods that were picked up when arriving back at the terminal (one can multiply vehicle capacity Qv with a factor between 0 and 1 to compensate for vehicles not being fully loaded on average). ? Then the total number of vehicle routes is equal to i∈O j∈D qij (ω)/Qv . Let the total region be denoted by A ? Suppose that each terminal serves origins and/or destinations in most of A (this is the case even ? with area ?(A). for distribution strategies such as one-terminal-nearest-terminal described in Section 2). The average number of vehicle routes per terminal is equal to i∈O j∈D qij (ω)/(N Qv ). Thus, if di?erent vehicle routes from the same terminal do not overlap, then the average area served per vehicle route is equal to ?(A) = ? ?(A)N Qv i∈O j∈D qij (ω)

Note that ?(A) depends on N and ω, but the notation does not indicate the dependence. Next we calculate the average number of delivery stops and the average number of pickup stops on a vehicle route, as a function of N , Q, and ω. The number of vehicle stops at an origin i ∈ O is at least max Ni , , and the number of vehicle stops at a destination j ∈ D is at least max Nj , j∈D qij (ω)/Qv i∈O qij (ω)/Qv (recall that Ni and Nj depend on N , Q, and ω). Thus, the average number of pickup stops per vehicle route is approximately i∈O max Ni , j∈D qij (ω)/Qv np = i∈O j∈D qij (ω)/Qv and the average number of delivery stops per vehicle route is approximately nd =
j∈D

max Nj ,
i∈O

i∈O qij (ω)/Qv

j∈D qij (ω)/Qv

Then, the approximate expected total detour distance D(N, Q) is calculated as follows: D(N, Q) :=
ω∈?

p(ω)

i∈O

j∈D qij (ω)

Qv

β

(np ? 1)?(A) +

(nd ? 1)?(A) +

?(A)

22

The approximation is obtained by substituting expressions for the average number of stops and the average area of the region served into the tour length approximation T n ≈ β n?(A). The reason np ? 1 and nd ? 1 are used in the tour length calculations is because the vehicle does not have to return to the ?rst pickup point or the ?rst delivery point after completing pickups or deliveries. The term β ?(A) approximates the average distance √ from the last delivery point to the ?rst pickup point on a vehicle route. Since x is a concave function, it follows from Jensen’s inequality that this overestimates the average tour length. In numerical experiments, this overestimation did not seem to have much of an e?ect on the selection of the optimal number of terminals. In practice, such overestimation is likely to be dominated by the amount by which actual route lengths on road networks di?er from the lengths according to simple metrics. Note that np is increasing in the numbers Ni of terminals serving origins i, nd is increasing in the numbers Nj of terminals serving destinations j, and D(N, Q) is increasing in np and nd . Also, the average area served per vehicle route ?(A) is increasing in the number N of terminals, and D(N, Q) is increasing in ?(A). Thus the expected total detour distance D(N, Q) is increasing in the number N of terminals and in the numbers of terminals serving origins and destinations, as it should be.

5.2.3

Search for Number of Terminals and Thresholds

Recall that we want to construct and solve an approximating problem (31):
N ∈{1,2,...}

min

? ? f (N ) := cN + V (N )

? taking into account the terminal ?xed cost cN and the approximate expected transportation cost V (N ). Without loss of generality, suppose that the unit of cost or the unit of distance has been scaled to make transportation ? cost per distance equal to 1. The approximate expected transportation cost V (N ) is given by minimizing the sum of the linehaul cost and the detour cost over the thresholds: ? V (N ) :=
0≤Q1 ≤···≤QN ?1

min

{L(N, Q) + D(N, Q)}

(36)

Note that the larger the value of N , the larger the set of thresholds that can be selected, and thus the smaller ? the value of V (N ). Note that if the total ?ows j∈D qij (ω) for all i ∈ O and i∈O qij (ω) for all j ∈ D are sorted, then all values of a threshold Qk between two successive sorted values of the total ?ow give the same values of L(N, Q) and D(N, Q). Thus, if N is small, say N ≤ 4, then problem (36) can easily be solved by enumerating all relevant values of the N ? 1 thresholds. If N is large, then problem (36) can be solved approximately by a neighborhood search on the set of relevant values of the N ? 1 thresholds. In addition, if a threshold Qk is changed from one interval in the sorted list of total ?ows to a neighboring interval, the resulting change in the values of L(N, Q) and D(N, Q) can be computed very quickly, because the value of Ni or Nj for only one origin i or destination j is a?ected by the change. ? Finally, problem minN ∈{1,2,...} cN + V (N ) can be solved by enumerating a range of reasonable values of N . ? The optimal value N is the number of terminals obtained with the CA method described above.

5.3

Terminal Location

Recall from Section 4 that it is not our purpose to model all the factors that should enter or do enter into the location of terminals. It is our purpose to test the CA method for selection of the number of terminals described in Section 5.2 by more detailed calculation of transportation costs. To facilitate such detailed calculation of vehicle routing costs, we need to calculate locations for the chosen number N of terminals. For that purpose,

23

we choose a set X of candidate locations, and solve the following problem: ? ? ? ? min cm um + (dim + dmj ) yij ?m u,? y ? ?
m∈X i∈O j∈D m∈X

subject to
m∈X yij ?m

yij ?m

= qij ?

for all i ∈ O, j ∈ D for all i ∈ O, j ∈ D, m ∈ X

≤ qij um ? um = N

m∈X yij ?m

≥ 0

for all i ∈ O, j ∈ D, m ∈ X for all m ∈ X

um ∈ {0, 1}

5.4

Choosing the Vehicle Fleet Sizes

The ?nal design decision required to test the CA method for selection of the number of terminals described in Section 5.2 is to choose the number of vehicles stationed at each terminal. Our approach is simple enumeration, using the operational decision procedures described in Section 6. First, for every scenario ω, we assign values to the decision variables zij (ω) that satisfy constraints (5)–(7), and such that the number of terminals that serve each origin i or destination j does not exceed the numbers Ni or Nj obtained from the thresholds Q and the ?ows q(ω) as described in Section 5.2. To do this, we use a simple heuristic described in Section 6.1. A lower bound on the number of vehicles required at each open terminal m is then given by Lm := max
ω∈? i∈O j∈D m zij (ω)qij (ω)

Qv

.

We use the detailed routing cost calculations described in Section 6.2 to calculate the cost of routing Lm , Lm + 1, Lm + 2, . . . , Lm + k vehicles from terminal m for some small number k. We then choose the number of vehicles for terminal m that minimizes the sum of vehicle cost and transportation cost from terminal m over the k +1 ?eet sizes. For each ?eet size, the detailed vehicle routing cost calculations for all scenarios described in Section 6.2 take a large amount of time. Therefore, as described in Section 5.2, the CA is based on an assumption of full vehicles, unlike the more detailed choice of vehicle ?eet sizes described in this section.

6

Operational Decisions

Section 5 described how to make the design decisions, namely selection of the number of terminals, location of the terminals, and choice of the vehicle ?eet sizes. In the process thresholds were also selected that control how many terminals serve an origin or destination, depending on the total ?ow from/to the origin/destination. In this section we describe methods for making operational decisions. Speci?cally, Section 6.1 describes a method for selecting which terminal to use for each origin-destination ?ow, and Section 6.2 describes a method for routing the vehicles from each terminal to do the pickups and deliveries.

6.1

Selection of Terminal for Each Origin-Destination Flow

This section describes how to decide through which terminal to route each origin-destination ?ow for a given set of ?ows q(ω). Ideally, for given open terminals u, vehicle ?eet sizes nv , and origin-destination ?ows q(ω), one would like to solve the integer linear program (4)–(7). However, solving (4)–(7) exactly does not seem to be practical. In this section, we describe a heuristic for choosing the decision variables z(ω) that specify through which terminal each origin-destination ?ow is routed. The heuristic does not require knowledge of the vehicle ?eet sizes nv , but does require knowledge of the thresholds Q.

24

Let Ni (k) and Nj (k) represent the sets of terminals serving origin i and destination j respectively at the kth iteration of the algorithm. We require that |Ni (k)| ≤ Ni , |Nj (k)| ≤ Nj , and Ni (k) ? U, Nj (k) ? U, where U := {m ∈ X : um = 1} denotes the set of open terminals. (Recall from Section 5.2 how the thresholds Q and the origin-destination ?ows q(ω) are used to determine the maximum numbers Ni and Nj of terminals serving origin i and destination j respectively.) For each origin i ∈ O, destination j ∈ D, and set N ? U, let m(i, j, N ) ∈ arg min {dim + dmj }
m∈N

denote a terminal from the set N that minimizes the distance from the origin i through a terminal in the set N to the destination j.

Algorithm to Select Terminal for Each Origin-Destination Flow:
(0) Initially, every origin and destination is served by only the center terminal, so that Ni (0) = Nj (0) := {0} for all i ∈ O and j ∈ D. (1) For each i ∈ O and j ∈ D, (a) If |Ni (k)| < Ni and |Nj (k)| < Nj , then set Ni+ (k) := Ni (k) ∪ m(i, j, U) (b) Else if |Ni (k)| < Ni and |Nj (k)| = Nj , then set Ni+ (k) := Ni (k) ∪ m(i, j, Nj (k)) (c) Else if |Ni (k)| = Ni and |Nj (k)| < Nj , set Ni+ (k) := Ni (k) and Nj+ (k) := Nj (k) ∪ m(i, j, Ni (k)); and Nj+ (k) := Nj (k); and Nj+ (k) := Nj (k) ∪ m(i, j, U);

(d) Else |Ni (k)| = Ni and |Nj (k)| = Nj . Set Ni+ (k) := Ni (k) (Note that |Ni+ (k)| ≤ Ni If and and Nj+ (k) := Nj (k).

|Nj+ (k)| ≤ Nj for all i ∈ O and j ∈ D.) and Nj+ (k) = Nj (k)

Ni+ (k) = Ni (k) (2) Choose (i , j )(k) ∈ arg max qij (ω)
i∈O,j∈D

for all i ∈ O and j ∈ D, then terminate the algorithm.

m∈Ni (k)∩Nj (k)

min

(dim + dmj ) ?

+ + m∈Ni (k)∩Nj (k)

min

(dim + dmj ) .

Set Ni (k + 1) := Ni+ (k) For all i ∈ O \ {i } and j ∈ D \ {j }, set Ni (k + 1) := Ni (k) Return to step 1. At the completion of the algorithm at ?nite iteration k , let
m zij (ω) = I{m=m(i,j,Ni (k )∩Nj (k ))} .

and Nj (k + 1) := Nj+ (k).

and Nj (k + 1) := Nj (k).

25

6.2

Vehicle Routing with Backhauls, Split Pickups and Deliveries

In this section we describe a method for routing the vehicles from each terminal to do given pickups and deliveries. The vehicle routes also provide a more accurate estimate of the transportation cost resulting from a given network design. The method described in this section is also, as far as we know, the ?rst heuristic proposed for the vehicle routing problem with backhauls and split pickups and deliveries, and thus may be of interest in itself. The vehicle routing problem with backhauls and split pickups and deliveries (VRPBS) was given in (8)–(18). Recall that the problem input is a terminal indexed by 0, a set O of origins, a set D of destinations, quantities to be picked up and delivered given by Q , vehicle movement costs between the origins, destinations, and the considered terminal given by d , and nv vehicles with capacity Qv each. Note that the origin-destination ?ows q(ω) and the output of the methods described in the previous sections provide the input for the vehicle routing problem, except for the ?eet size nv , which can be determined by enumeration, as described in Section 5.4. We describe a cluster-?rst, route-second heuristic for the VRPBS that uses various ideas of the heuristic for the vehicle routing problem with backhauls (VRPB) proposed by Toth and Vigo (1999). In both the VRPB and the VRPBS, deliveries have to be performed before pickups on the same route. The heuristic of Toth and Vigo (1999) has to be modi?ed for the following reasons. First, as already pointed out, in the VRPBS multiple vehicles are allowed to visit each origin and destination, whereas in the VRPB exactly one vehicle must visit each origin and destination. One reason this modi?cation is needed is because it often holds that Qi > Qv for some origins or destinations i ∈ O ∪ D , and thus some origins and destinations must be visited by more than one vehicle (recall the typical nonuniform distribution of pickup and delivery quantities shown in Figure 6). Second, in the VRPB considered by Toth and Vigo (1999), each of the nv vehicles must visit at least one destination, and thus no vehicle may travel directly from the terminal to an origin. In the VRPBS, fewer than nv vehicles may visit origins or destinations, and vehicles may travel directly to origins (and thus not visit any destinations).

6.2.1

Initial Splitting Step

As pointed out, for some points i ∈ O ∪ D , it may hold that Qi > Qv . If point i has quantity Qi > Qv , then it should be served by multiple vehicles, and all these vehicles except possibly one should carry a full load associated with this point, and the remaining quantity should be carried by another vehicle which can combine the remaining quantity with additional loads associated with other points. For example, if a destination j has Qj = 3.5Qv , then each of 3 vehicles should deliver full vehicle loads to j and afterward move to pick-up points, and one vehicle should deliver a load of size 0.5Qv to j and the rest of the vehicle’s space could be used for other deliveries on its route. Thus, in the ?rst step each point i and its quantity Qi is split, after which each new point i has quantity Qi ≤ Qv . Speci?cally, for each point i ∈ O ∪ D , create Qi /Qv copies of point i, of which Qi /Qv new points i have quantity Qi = Qv , and if Qi /Qv < Qi /Qv , then the remaining point i has quantity Qi = Qi ? Qi /Qv Qv . Let O denote the set of new origins, D denote the set of new destinations, V := {0} ∪ O ∪ D denote the new set of nodes, A := (i, j) ∈ (V )2 \ O × D : i = j denote the new set of arcs, and Qi and Qj denote the new quantities to be picked up and delivered. For (i, j) ∈ A , the vehicle movement costs di,j are obtained from the given vehicle movement costs d in the obvious way, except if i and j correspond to the same original point, in which case di,j := ε for some chosen ε > 0. Next one can de?ne the following VRPB with input data O , D , V , A , Q , d , nv , and Qv . The decision variables are xij nO v nD v := := := 1 0 if a vehicle travels on arc (i, j) otherwise

number of vehicles performing pickups number of vehicles performing deliveries

Note that the number of vehicles used is given by max{nO , nD }. For each set S ? O or S ? D , let σ(S) v v denote the minimum number of vehicles needed to serve each point in S if exactly one vehicle must visit each

26

point, that is, σ(S) is the optimal value of the bin packing problem with item sizes given by Qi for i ∈ S and bin size Qv . For each combination of nO ∈ v let τ O , D , Q , d , nO , nD v v := minx
(i,j)∈A i∈O

Qi /Qv , . . . , nv dij xij xji = 1

and nD ∈ v

j∈D

Qj /Qv , . . . , nv , (37)

subject to
{j : (j,i)∈A }

for all i ∈ O ∪ D for all i ∈ O ∪ D

(38) (39) (40) (41) (42) (43) (44) (45) (46)

xij
{j : (i,j)∈A }

= 1

xi0 = nO v
i∈O

x0i = max{0, nO ? nD } v v
i∈O

x0j
j∈D

= nD v

xj0 = max{0, nD ? nO } v v
j∈D

xij
{(i,j)∈A : i∈S,j ∈S} /

≥ σ(S) for all S ? O , |S| ≥ 2 ≥ σ(S) for all S ? D , |S| ≥ 2 for all (i, j) ∈ A

xij
{(i,j)∈A : i∈S,j∈S} /

xij

∈ {0, 1}

denote the optimal value of the VRPB given that exactly nO vehicles serve the origins O and exactly nD vehicles v v serve the destinations D . If nO or nD is too small so that (38)–(46) is infeasible, then τ O , D , Q , d , nO , nD := v v v v ∞. Then the VRPB is given by τ (O , D , Q , d , nv ) := minnO ,nD v v subject to cv max{nO , nD } + τ v v nO v nD v ∈
i∈O

O , D , Q , d , nO , nD v v

Qi /Qv , . . . , nv

?? ? ? ? ? ? ∈ Qj /Qv ? , . . . , nv ? ?? ? ? ?j∈D

6.2.2

Initial Clustering Step

Note that, for any i ∈ O , if Qi = Qv , then (38)–(46) imply that i is visited by exactly one vehicle that does not visit any other origin. Similarly, for any j ∈ D , if Qj = Qv , then j is visited by exactly one vehicle that ? ? does not visit any other destination. Let O := {i ∈ O : Qi = Qv } and D := {j ∈ D : Qj = Qv }. Then for ? , (38) can be replaced by all i ∈ O x0i +
j∈D

xji = 1 for all j ∈ O

(47) (48)

xji = 0 ? and for all j ∈ D , (38) can be replaced by x0j xij = 1 = 0

(49) for all i ∈ D (50)

27

? Also, for all i ∈ O , (39) can be replaced by xi0 = 1 xij ? and for all j ∈ D , (39) can be replaced by xj0 +
i∈O

(51) for all j ∈ O (52)

= 0

xji = 1 for all i ∈ D
? i∈O \O

(53) (54) ? ? xi0 = nO ? |O |. If O \ O v = ?,

xji = 0

Note that, given (51), constraint (40) holds if and only if O ? then it follows from nO ≥ v i∈O Qi /Qv that nv ? |O | ≥ nD v

constraint (42) holds if and only if ? j∈D \D x0j = D ? |≥ nv ? |D ? j∈D \D Qj /Qv ≥ 1. Also note that constraints (47) and (51) imply that constraint (44) holds ? . Similarly, constraints (49) and (53) imply that constraint (45) holds for all S ? D . Furthermore, ? for all S ? O ? , and not for S ? O such constraints (48) and (52) imply that constraint (44) is required only for S ? O \ O ? ? that S ∩ O = ? and S ∩ (O \ O ) = ?. Similarly, constraints (50) and (54) imply that constraint (45) is ? , and not for S ? D such that S ∩ D = ? and S ∩ (D \ D ) = ?. ? ? required only for S ? D \ D Note that it follows from constraints (38)–(43) that j∈D i∈O xji = j∈D i∈D xji ? / j∈D xj0 = xij ? j∈D xj0 = j∈D x0j ? j∈D xj0 = nD ? max{0, nD ? nO } = min{nO , nD }. Thus we v v v v v i∈D / j∈D can add the following redundant constraint: xji = min{nO , nD } v v
j∈D i∈O

Qi /Qv ≥ 1. Similarly, given (49), ? ? ? |D |, and if D \ D = ?, then it follows that
? i∈O \O

(55)

Next, we add the following redundant constraints that follow from constraints (38) and (39) respectively: x0i +
j∈D

xji xji
i∈O

≤ 1 ≤ 1

? for all i ∈ O \ O ? for all j ∈ D \ D

(56) (57)

xj0 +

Also, we add the following redundant constraints that follow from constraints (44) and (45) respectively: xij
{(i,j)∈A : i∈S,j ∈S} /

≥ 1 ≥ 1

? for all S ? O \ O , |S| ≥ 2 ? for all S ? D \ D , |S| ≥ 2

(58) (59)

xij
{(i,j)∈A : i∈S,j∈S} /

Next we formulate a Lagrangian relaxation for problem (37)–(59). Let the multipliers associated with con? ? straint (38) for i ∈ O \ O and associated with constraint (39) for i ∈ D \ D be denoted by λi , and let the ? ? multipliers associated with constraint (44) for S ? O \ O and associated with constraint (45) for S ? D \ D

28

be denoted by ?S ≥ 0. Then the corresponding Lagrangian relaxation is as follows: ? L O , D , Q , d , nO , nD , λ, ? v v := minx
(i,j)∈A

? xji ? 1?

dij xij + ? λi ?
? i∈D \D ? i∈O \O

λi ?
{j : (j,i)∈A }

?

+

xij ? 1?
{j : (i,j)∈A }

?

? xij ?
: i∈S,j ∈S} /

+
? S?O \O

?S ?σ(S) ? ? ?S ?σ(S) ?
? S?D \D {(i,j)∈A : i∈S,j∈S} / {(i,j)∈A

? xij ? (60) (61) (62) (63) (64) (65) (66) (67) (68) (69) (70) (71) (72) (73) (74) (75) (76) (77)

+ subject to x0i +
j∈D

xji = 1 xji ≤ 1
j∈D

? for all i ∈ O ? for all i ∈ O \ O ? for all i ∈ O , j ∈ O ? for all j ∈ D

x0i + xij x0j xij = 0 = 1 = 0

xij
{i : (i,j)∈A }

= 1

? for all i ∈ D , j ∈ D ? for all j ∈ D \ D ? for all i ∈ O ? for all i ∈ O , j ∈ O ? for all i ∈ O \ O ? for all j ∈ D ? for all j ∈ D \ D ? for all i ∈ D , j ∈ D

xi0 = 1 xij = 0 xij
{j : (i,j)∈A }

= 1

xj0 +
i∈O

xji = 1 xji ≤ 1
i∈O

xj0 + xij = 0

? xi0 = nO ? O v
? i∈O \O

x0i = max{0, nO ? nD } v v
i∈O

x0j
? j∈D \D

? = nD ? D v

xj0 = max{0, nD ? nO } v v
j∈D

xji = min{nO , nD } v v
j∈D i∈O

xij
{(i,j)∈A : i∈S,j ∈S} /

≥ 1 ≥ 1

? for all S ? O \ O , |S| ≥ 2(78) ? for all S ? D \ D , |S| ≥ 2(79) for all (i, j) ∈ A (80)

xij
{(i,j)∈A : i∈S,j∈S} /

xij

∈ {0, 1}

29

For each (i, j) ∈ A , let ? ? dij ? ? ? d ? ij ? ? ? ? dij ? ? ? ? d ij ? := dij ? dij ? ? ? d ? ij ? ? ? ? d ? ij ? ? ? dij

+ λj ? {S?D \D : j∈S} ?S ? ? {S?O \O : i∈S} ?S ? + λi + λj ? {S?O \O : i∈S,j ∈S} ?S ? / + λi ? {S?D \D : i∈S,j∈S} ?S ? / + λi + λj

if if if if if if if if

? ? i, j ∈ {0} ∪ O ∪ D ? ? ? i ∈ {0} ∪ O ∪ D , j ∈ O \ O ? ? ? i ∈ {0} ∪ O ∪ D , j ∈ D \ D ? ? i ∈ O \ O , j ∈ {0} ∪ O ? ? ? i ∈ D \ D , j ∈ {0} ∪ O ∪ D ? i, j ∈ O \ O ? i, j ∈ D \ D ? ,j ∈ O \ O ? i∈D \D

(81)

Then the objective function (60) is equal to ? dij xij ?
(i,j)∈A ? i∈O \O

λi ?
? i∈D \D

λi +
? S?O \O

?S σ(S) +
? S?D \D

?S σ(S)

To interpret feasible solutions of the Lagrangian relaxation (60)–(80), partition the arcs A into the following subsets: A1 A2 A3 A4 A5 := := := := ? (i, 0) : i ∈ O ? ∪ (0, j) : j ∈ D ? ∪ (i, j) : i ∈ O , j ∈ O ? ∪ (i, j) : i ∈ D , j ∈ D ? ∪ (i, j) : i ∈ D , j ∈ D

? (i, j) : i ∈ O , j ∈ O

? ? (i, j) : i ∈ O \ O , j ∈ {0} ∪ O \ O ? ? (i, j) : i ∈ {0} ∪ D \ D , j ∈ D \ D

:= {(i, j) : i ∈ {0} ∪ D , j ∈ {0} ∪ O }

First, observe that each of the constraints (61)–(79) involves decision variables xij for arcs (i, j) in only one of the subsets above. Speci?cally, constraints (64) and (67) involve arcs in A1 only; constraints (63), (65), (68), and (72) involve arcs in A2 only; constraints (69), (73), and (78) involve arcs in A3 only; constraints (66), (75), and (79) involve arcs in A4 only; and constraints (61), (62), (70), (71), (74), (76), and (77) involve arcs in A5 only. Also, clearly each individual constraint in (80) involves an arc in only one of the subsets above. Thus, the Lagrangian relaxation (60)–(80) decomposes into 5 subproblems, corresponding to the sets of arcs and constraints identi?ed above. Next we consider each of these 5 subproblems in turn. For the subproblem involving A1 , it follows from constraints (64) and (67) that xij = 1 for all (i, j) ∈ A1 . For the subproblem involving A2 , it follows from constraints (63), (65), (68), and (72) that xij = 0 for all (i, j) ∈ A2 . For the subproblem involving A3 , it follows from constraint (69) that exactly one arc out of each node ? ? ? i ∈ O \ O must be chosen, it follows from constraint (73) that exactly nO ? O arcs from nodes i ∈ O \ O v ? to the terminal node 0 must be chosen, and it follows from constraint (78) that for each subset S ? O \ O , at least one arc out of the subset must be chosen, and thus it follows from constraints (69), (73), and (78) ? that the chosen arcs may not form any cycles in O \ O . In other words, the chosen arcs in A3 must form a ? spanning anti-arborescence on the nodes {0} ∪ O \ O with the terminal node 0 being the root node and with O ? | at node 0. Thus the subproblem involving A is a shortest spanning anti-arborescence ?xed indegree nv ? |O 3 ? problem with ?xed indegree K O := nO ? |O | at the terminal node 0 (K O -SSAA), and with arc costs given by v ? dij . Similarly, it follows from constraints (66), (75), and (79) that the subproblem involving A4 is a shortest ? spanning arborescence problem with ?xed outdegree K D := nD ? |D | at the terminal node 0 (K D -SSA), and v ? . Problems K O -SSAA and K D -SSA can be solved in O(|O \ O |2 ) and O(|D \ D |2 ) ? ? with arc costs given by dij time respectively, for example with the algorithm of Gabow and Tarjan (1984).

30

Next we show that the subproblem involving A5 can be represented as a network ?ow problem. The network ?ow problem has a node for each node in V , as well as an additional source node s and sink node t. The supply at the source node and the demand at the sink node are both equal to max{nO , nD }. There is an arc (s, j) from v v ? the source node s to each node j ∈ D with cost 0. The lower bound of the ?ow on each arc (s, j) is 1 if j ∈ D ? and 0 if j ∈ D \ D . The upper bound of the ?ow on each arc (s, j) is 1. There is an arc (j, i) from each node ? j ∈ D to each node i ∈ O with cost dji . The lower bound of the ?ow on each such arc (j, i) is 0, and the upper bound of the ?ow on each such arc (j, i) is 1. There is an arc (i, t) from each node i ∈ O to the sink ? ? node t with cost 0. The lower bound of the ?ow on each arc (i, t) is 1 if i ∈ O and 0 if i ∈ O \ O , and the O D upper bound of the ?ow on each arc (i, t) is 1. Suppose that nv > nv . Then there is an arc (s, 0) with cost 0. The lower bound and the upper bound of the ?ow on arc (s, 0) are both equal to nO ? nD . There is also an v v ? arc (0, i) from node 0 to each node i ∈ O with cost d0i . The lower bound of the ?ow on each such arc (0, i) is 0, and the upper bound of the ?ow on each such arc (0, i) is 1. If nO < nD , then there is an arc (j, 0) from v v ? each node j ∈ D to node 0 with cost dj0 , lower bound 0, and upper bound 1, and an arc (0, t) with cost 0, and lower bound and upper bound both equal to nD ? nO . It is easy to see that for every solution that satis?es v v constraints (61), (62), (70), (71), (74), (76), (77), and (80) (for (i, j) ∈ A5 ), there is a feasible integer ?ow for the network ?ow problem described above with the same cost, and vice versa. The network ?ow problem can be solved in O(max{nO , nD }(|O | + |D |)2 ) time, for example with a shortest augmenting path algorithm; see, for v v example, Ahuja et al. (1993). ? ? Next we brie?y address the following two issues. First, the number of subsets S ? O \ O and S ? D \ D in the Lagrangian objective (60) may be very large. Second, we want to ?nd multipliers λ, ? that solve the Lagrangian dual problem max L O , D , Q , d , nO , nD , λ, ? : ? ≥ 0 v v
λ,?

We use the same subgradient optimization procedure described in Toth and Vigo (1997) to simultaneously ? ? address both issues. Brie?y, instead of enumerating all subsets S ? O \ O and S ? D \ D , at each O D major iteration the procedure identi?es subtrees in the K -SSAA and K -SSA that violate constraints (44) and (45) respectively, and adds the terms corresponding to the subsets S of nodes in the violating subtrees to ? the Lagrangian objective (60). The identi?cation of the violating subtrees can be done in O(|O \ O |) and ? |) time respectively. During each major iteration, the multipliers λ and ? are updated by performing O(|D \ D several minor iterations of a subgradient search. For additional details, we refer to Toth and Vigo (1997).

6.2.3

Second Splitting Step

When progress made by the subgradient optimization procedure in solving the Lagrangian dual problem slows down, the procedure is stopped. The last K O -SSAA and K D -SSA constructed are considered. If there are no subarborescences in the K O -SSAA and K D -SSA that violate constraints (44) and (45) respectively, then set ? ? ? ? ? ? ? ? O := O , O \ O := O \ O , D := D , D \ D := D \ D , and continue with the assignment-routing step described in Section 6.2.5. Otherwise, a second round of node splitting is performed, as described in this section. ? Suppose there are subarborescences in the K D -SSA constructed on nodes {0} ∪ D \ D that violate con? is partitioned into K D := nD ? |D | subsets D , k = 1, . . . , K D , corresponding ? straint (45). The set D \ D v k ? to the K D subarborescences in the K D -SSA rooted at node 0. For each j ∈ D \ D and k ∈ {1, . . . , K D }, let djk := min{dji : i ∈ Dk } denote the distance between node j and subarborescence k. Note that if j ∈ Dk , then djk = 0. Next the following transportation problem is solved to assign loads to subarborescences in such a way that the total load assigned to each subarborescence is less than the vehicle capacity, and such that the total

31

assignment cost is minimized.
KD

miny
? j∈D \D k=1 KD

djk yjk

(82)

subject to
k=1

yjk

= Qj yjk ≤ Qv

? for all j ∈ D \ D for all k ∈ {1, . . . , K D } ? for all j ∈ D \ D , k ∈ {1, . . . , K D }

(83) (84) (85)

? j∈D \D

yjk

≥ 0

Let y ? denote an optimal solution of problem (82)–(85). The new set D of nodes and their loads Qi , i ∈ D ? ? ? ? are determined as follows: D := D and Qi := Qv for all i ∈ D . Also, for each j ∈ D \ D and each D ? ? with load size Q := y ? . If there are k ∈ {1, . . . , K } such that yjk > 0, there is a node i ∈ D \ D i jk ? ? ? ? no subarborescences in the K D -SSA that violate constraint (45), then D := D , D \ D := D \ D , and ? , O \ O , and the loads Q for i ∈ O , are determined similarly ? Qi := Qi for all i ∈ D . The sets O i based on the K O -SSAA. Note that after this second splitting, we are guaranteed to ?nd a feasible solution of D the VRPBS, as long as nO ≥ v i∈O Qi /Qv and nv ≥ j∈D Qj /Qv . Let V := {0} ∪ O ∪ D denote the new set of nodes, and A := (i, j) ∈ (V )2 \ O × D : i = j denote the new set of arcs. For (i, j) ∈ A , the vehicle movement costs di,j are obtained from the given vehicle movement costs d, and as before if i and j correspond to the same original point, then di,j := ε for some chosen ε > 0.

6.2.4

Second Clustering Step

After the second splitting step, a second clustering step is performed. The second clustering step is the same as the initial clustering step described in Section 6.2.2, but with O , D , V , A , Q , and d instead of O , D , V , A , Q , and d respectively.

6.2.5

Assignment-Routing Step

Consider the K O -SSAA and K D -SSA constructed in the last clustering step. The set D is partitioned into ? ? nD subsets as follows. First, D is partitioned into its singletons: for each k = 1, . . . , |D |, choose without v ? and set D := {i}. Second, the set D \ D is partitioned into K D := nD ? |D | ? ? replacement i ∈ D v k D D D ? subsets Dk , k = |D | + 1, . . . , nv , corresponding to the K subarborescences in the K -SSA rooted at node 0. Similarly, set O is partitioned into nO subsets Ol , l = 1, . . . , nO . For each k ∈ {1, . . . , nD } and v v v l ∈ {1, . . . , nO }, consider the traveling salesman problem (TSP) with node set Vkl := {0} ∪ Dk ∪ Ol , arc v set Akl := (i, j) ∈ (Vkl )2 \ Ol × Dk : i = j , and arc costs dij , (i, j) ∈ Akl . Note that by construction the TSP has the precedence constraint that after 0, all nodes in Dk must be visited before any nodes in Ol are visited, and is thus called the traveling salesman problem with backhauls (TSPB). Let dkl denote the optimal objective value, or an estimate of the optimal objective value, of the TSPB. Similar to Toth and Vigo (1999), we use the farthest insertion heuristic to obtain a solution of the TSPB, and set dkl equal to the objective value of the solution produced by the heuristic. In addition, for each k ∈ {1, . . . , nD }, consider the TSP with v node set Vk0 := {0} ∪ Dk , arc set Ak0 := (i, j) ∈ (Vk0 )2 : i = j , and arc costs dij , (i, j) ∈ Ak0 . Let dk0 denote the optimal objective value, or an estimate of the optimal objective value, of the TSP. Similarly, for each l ∈ {1, . . . , nO }, consider the TSP with node set V0l := {0} ∪ Ol , arc set A0l := (i, j) ∈ (V0l )2 : i = j , and v arc costs dij , (i, j) ∈ A0l . Let d0l denote the optimal objective value, or an estimate of the optimal objective value, of the TSP.

32

Next the following assignment problem is solved to combine subsets of D
nD v nO v nD nO v v

and O

into vehicle routes. (86)

minz
k=1

dk0 zk0 +
l=1 nO v

d0l z0l +
k=1 l=1

dkl zkl for all k ∈ {1, . . . , nD } v for all l ∈ {1, . . . , nO } v

subject to

zk0 +
l=1 nD v

zkl = 1

(87)

z0l +
k=1

zkl = 1

(88)

zk0 , z0l , zkl ∈ {0, 1}

for all k ∈ {1, . . . , nD }, l ∈ {1, . . . , nO } (89) v v

? Let z ? denote an optimal solution of problem (86)–(89). If zkl = 1, then the nodes in {0}, Dk and Ol are ? included in a vehicle route given by the TSPB discussed above. Similarly, if zk0 = 1, then the nodes in {0} ? and Dk are included in a vehicle route, and if z0l = 1, then the nodes in {0} and Ol are included in a vehicle route, given by the TSP. The resulting routes may violate vehicle capacity constraints. The improvement step discussed in the next section attempts to modify routes to eliminate violations of capacity constraints.

6.2.6

Improvement Heuristic

After construction of the vehicle routes as described in the previous section, an improvement heuristic is applied to modify routes to eliminate violations of capacity constraints, and to reduce the costs. (Note that another reasonable option is to ?rst obtain a feasible solution by applying the assignment-routing step in Section 6.2.5 to the feasible solution obtained in the second splitting step in Section 6.2.3, and to then use an improvement heuristic to ?nd a solution with a better objective value.) The improvement heuristic is the same as the postoptimization procedure described in Toth and Vigo (1999). Brie?y, intra-route 2-exchanges and 3-exchanges are performed in each route to reduce the cost of the route. Also, inter-route 1-exchanges and 2-exchanges are performed to reduce the amount by which capacity constraints are violated, and to reduce the costs of the routes. We refer to Toth and Vigo (1999) for additional details. Let τ O , D , Q , d , nO , nD denote the total cost of the resulting solution of the VRPBS, given that exactly ? v v nO vehicles serve the origins O and exactly nD vehicles serve the destinations D . (Note that O , D , V , v v A , Q , and d are determined by O , D , Q , and d , and thus we may denote τ as a function of O , D , Q , ? and d .) Then let τ (O , D , Q , d , nv ) := ? minnO ,nD v v subject to cv max{nO , nD } + τ O , D , Q , d , nO , nD ? v v v v nO v nD v ∈
i∈O

Qi /Qv , . . . , nv

?? ? ? ? ? ? ∈ Qj /Qv ? , . . . , nv ? ?? ? ? ?j∈D

denote the total cost of the resulting solution of the VRPBS. Recall that, as referred to in Section 5.4, the number nm of vehicles for each terminal m can be chosen by v solving minnm v subject to Cv nm + v
ω∈?

p(ω)? (Om (z m (ω)), Dm (z m (ω)), Qm (z m (ω), ω), dm (z m (ω)), nm ) τ v

nm ∈ {Lm , Lm + 1, Lm + 2, . . . , Lm + k} v

where z m (ω) can be chosen as described in Section 6.1, and Om (z m (ω)), Dm (z m (ω)), Qm (z m (ω), ω), and dm (z m (ω)) are de?ned in Section 3.

33

7

Evaluation of the Proposed CA Solutions

To evaluate our solutions, we considered eight combinations of small and large datasets, with origins and destinations uniformly and nonuniformly distributed, and with origin-destination ?ows correlated and uncorrelated, as described in more detail below. For each set of origin and destination locations, we generated ten ?ow scenarios according to a speci?ed distribution as described below, and used these ten design scenarios to make the design decisions, namely selection of the number of terminals and the thresholds, and the location of the terminals, using the method described in Section 5. We also solved the TEMC problem (24)–(28) using the same ten design scenarios, which provided a design for comparison. We selected the same vehicle ?eet sizes, namely Lm + 2, to evaluate both the CA solution and the TEMC solution. Using the same ?ow distribution, we then generated ?ow scenarios for evaluation, independent of the design scenarios, and calculated the cost for each evaluation scenario using the approach described in Section 6 to make the operational decisions. Speci?cally, for each evaluation scenario we selected which terminal served each origin-destination ?ow and how to route the vehicles from each terminal to do the pickups and deliveries. Because the time involved in performing the detailed routing calculations described in Section 6.2 was large, it was not practical to generate a huge number of evaluation scenarios and calculate the cost associated with each one. Therefore, we only generated enough scenarios (around 35) to ensure that the sample standard deviation of the sample average cost di?erence was less than half the sample average cost di?erence. Speci?cally, let {ω1 , . . . , ωM } cCA (ωm ) cT EM C (ωm ) := := := set of evaluation scenarios cost of the CA design for scenario ωm cost of the TEMC design for scenario ωm .

Then, the evaluation sample size M is chosen large enough so that ? 1 1 ? 2 (cCA (ωm ) ? cT EM C (ωm )) ? ? M (M ? 1) m=1 < 1 2
M m=1 M M m=1

(cCA (ωm ) ? cT EM C (ωm )) M

2

? ? ?

(cCA (ωm ) ? cT EM C (ωm )) . M

(90)

The origin and destination locations. The industry dataset with 148 origins and 36 destinations, shown
in Figure 5, provided a set of origins and destinations for evaluation, namely for the two cases (correlated and uncorrelated ?ows) with small datasets and nonuniform distributions of origins and destinations. To generate origin and destination locations “uniformly”, we selected 1097 cities roughly uniformly distributed across the USA. For the two cases with small datasets and uniform distributions of origins and destinations, we sampled 148 of these 1097 cities without replacement to be origins, and 41 of these 1097 cities without replacement to be destinations. For the two cases with large datasets and uniform distributions of origins and destinations, the numbers of origins and destinations chosen were 243 and 105 respectively. Finally, for the two cases with large datasets and nonuniform distributions of origins and destinations, we used only the 788 of the 1097 cities that were approximately west of Colorado and east of Illinois, and then chose 248 origins and 134 destinations from the 788 cities without replacement. The reason for having more origins than destinations was to mimic the observed situations in the motivating industrial datasets.

Flow generation. Next we describe how we generated ?ows between origins and destinations for each of
the eight cases. We ?rst decided which origin-destination pairs would have positive ?ow. Because long distance positive ?ows tend to be less frequent than positive ?ows between origin-destination pairs located closer together, we selected each origin-destination pair with the origin located approximately west of Colorado (east of Illinois) and with the destination located approximately east of Illinois (west of Colorado) with probability 0.15 to have

34

Uncorrelated Demands Origins and Destinations Total number Location Distribution Small Uniform Nonuniform (Industry dataset) Large Uniform Nonuniform Correlated Demands Origins and Destinations Total number Location Distribution Small Uniform Nonuniform (Industry dataset) Large Uniform Nonuniform

Number of Terminals TEMC CA 2 1 3 1 2 3 3 4

Cost TEMC 644,489 382,792 873,348 1,018,598

CA 631,008 353,982 858,519 990,303

STD 5150 5688 5288 8890

% Decrease 2.09% 7.53% 1.70% 2.78%

Number of Terminals TEMC CA 2 1 3 1 3 1 3 4

Cost TEMC 688,519 616,531 675,098 989,678

CA 673,474 589,997 599,387 965,099

STD 5124 5333 5383 9691

% Decrease 2.19% 4.30% 11.21% 2.48%

Table 2: Numbers of terminals chosen and resulting costs.
zero ?ow, independently for all origin-destination pairs. In addition, for the four large cases, each remaining origin-destination pair, including the above mentioned “cross-country” pairs as well as pairs closer to each other, was assigned zero ?ow with probability 0.5, again independently for all origin-destination pairs. The remaining origin-destination pairs were assigned positive ?ows as described below. The expected number of positive ?ows was thus 1 I(“large” case) × 2 total number of origin-destination pairs ?0.15 × the number of “cross-country” origin-destination pairs ,

which yielded approximately 5000 positive ?ows for the small cases and 10,000 positive ?ows for the large cases. For the origin-destination pairs having positive ?ow, we generated the size of the ?ow according to a uniform distribution between 0 and an origin-speci?c upper bound that follows approximately an exponential distribution among origins. Recall from Figure 6 that the industrial dataset shows that a small fraction of the origindestination pairs accounts for a large fraction of the total ?ow. Similarly, a small fraction of the origins, mostly the origins of bulky goods, accounts for a large fraction of the total ?ow. To mimic this skewed ?ow distribution, we generated the ?ow upper bounds of the origins as follows. We ?rst generated a random permutation of the origins to label the origins from 1 to |O|. Then the ?ow from each origin i to any destination with positive origin-destination ?ow was distributed uniformly on [0, 1500 exp (??(i ? 1))] , for i ∈ {1, . . . , |O|}. Here, ? was chosen as 0.06 for the small uniform cases, 0.055 for the small nonuniform cases, 0.045 for the large uniform cases, and 0.043 for the large nonuniform cases. Finally, for the four cases with uncorrelated ?ows, we generated independent random variables with uniform distributions as described above for every origin-destination pair with positive ?ow. For the four cases with correlated ?ows, we ?rst generated independent uniform random variables as described above, and then multiplied all these independently generated ?ows by a single realization of a uniform (0.7, 1.3) random variable. It was of interest to also test our approach with correlated ?ow distributions because in practice ?ows are often correlated. For example, for home improvement retailers, shipments between origins and destinations tend to be relatively large during spring and summer.

Results. Table 2 displays the number of terminals chosen in both the CA and TEMC solutions, the average
cost of both solutions, the sample standard deviation of the sample average cost di?erence given in (90), and

35

the percentage cost improvement of the CA solution over the TEMC solution. All the cases investigated showed cost improvement of at least 1%, with one case showing a cost improvement of more than 10%. Observe that the cases with the largest percentage improvements (4.30%, 7.53%, and 11.21%) are also the cases in which the number of terminals chosen by the CA vs. the TEMC solution di?ers by more than 1.

Acknowledgements
This research was supported in part by NSF grant DMI-0113881. We also thank Pete Viehweg for his support of this work, and Brynn Conover for her help with the data.

References
R. K. Ahuja, T. L. Magnanti, and J. B. Orlin. Network Flows. Prentice-Hall, Englewood Cli?s, NJ, 1993. S. Anily. The vehicle-routing problem with delivery and back-haul options. Naval Research Logistics Quarterly, 43:415–434, 1996. C. Archetti, M. G. Speranza, and A. Hertz. A tabu search algorithm for the split delivery vehicle routing problem. Transportation Science, 40(1):64–73, 2006. J. Beardwood, J. H. Halton, and J. M. Hammersley. The shortest path through many points. Proceedings of the Cambridge Philosophical Society, 55:299–327, 1959. M. Beckmann. A continuous model of transportation. Econometrica, 20(4):643–660, 1952. D. E. Blumenfeld, L. D. Burns, C. F. Daganzo, M. C. Frick, and R. W. Hall. Reducing logistics costs at general motors. Interfaces, 17(1):26–47, 1987. D. E. Blumenfeld, L. D. Burns, J. D. Diltz, and C. F. Daganzo. Analyzing trade-o?s between transportation, inventory and production costs on freight networks. Transportation Research, 19B(5):361–380, 1985. L. D. Burns, R. W. Hall, D. E. Blumenfeld, and C. F. Daganzo. Distribution strategies that minimize transportation and inventory costs. Operations Research, 33(3):469–490, 1985. J. F. Campbell. Designing logistics systems by analyzing transportation, inventory and terminal cost tradeo?s. Journal of Business Logistics, 11(1):159–179, 1990a. J. F. Campbell. Freight consolidation and routing with transportation economies of scale. Transportation Research, 24B:345–361, 1990b. J. F. Campbell. Continuous and discrete demand hub location problems. Transportation Research, 27B(6): 473–482, 1993a. J. F. Campbell. One-to-many distribution with transshipments: An analytic model. Transportation Science, 27: 330–340, 1993b. C. F. Daganzo. An approximate analytic model of many-to-many demand responsive transportation systems. Transportation Research, 12B:325–333, 1978. C. F. Daganzo. The distance traveled to visit N points with a maximum of C stops per vehicle: An analytic model and an application. Transportation Science, 18:331–350, 1984a. C. F. Daganzo. The length of tours in zones of di?erent shapes. Transportation Research, 18B(2):135–145, 1984b.

36

C. F. Daganzo. Supplying a single location from heterogeneous sources. Transportation Research, 19B(5): 409–419, 1985. C. F. Daganzo. The break-bulk role of terminals in many-to-many logistic networks. Operations Research, 35: 543–555, 1987. C. F. Daganzo. A comparison of in-vehicle and out-of-vehicle freight consolidation strategies. Transportation Research, 22B:173–180, 1988a. C. F. Daganzo. Shipment composition enhancement at a consolidation center. Transportation Research, 22B: 103–124, 1988b. C. F. Daganzo and R. Hall. A routing model for pickups and deliveries: No capacity restrictions on the secondary items. Transportation Science, 27:315–329, 1993. C. F. Daganzo and G. F. Newell. Con?guration of physical distribution networks. Networks, 16:113–132, 1986. A. Dasci and G. Laporte. A continuous model for multistore competitive location. Operations Research, 53(2): 263–280, 2005. A. Dasci and V. Verter. A continuous model for production-distribution system design. European Journal of Operational Research, 129(2):287–298, 2001. M. Dror, G. Laporte, and P. Trudeau. Vehicle routing with split deliveries. Discrete Applied Mathematics, 50: 239–254, 1994. M. Dror and P. Trudeau. Savings by split delivery routing. Transportation Science, 23(2):141–145, 1989. H. N. Gabow and R. E. Tarjan. E?cient algorithms for a family of matroid intersection problems. Journal of Algorithms, 5:80–131, 1984. G. Ghiani, G. Laporte, and R. Musmanno. Introduction to Logistics Systems Planning and Control. John Wiley & Sons, Ltd., ChiChester, England, 2004. M. Goetschalckx and C. Jacobs-Blecha. The vehicle routing problem with backhauls. European Journal of Operational Research, 42(1):39–51, 1989. R. W. Hall. Travel distance through consolidation terminals on a rectangular grid. Journal of the Operational Research Society, 35:1067–1078, 1984. R. W. Hall. Comparison of strategies for routing shipments through transportation terminals. Transportation Research, 21A:421–429, 1987. R. W. Hall. Characteristics of multi-stop/multi-terminal delivery routes, with backhauls and unique items. Transportation Research, 25B:391–403, 1991. R. W. Hall. Design for local area freight networks. Transportation Research, 27B:79–95, 1993. L. Kantorovitch. On the translocation of masses. Comptes Rendus (Doklady) de l’Acad?mie des Sciences de e l’URSS, 37(7–8):199–201, 1942. reprinted in Management Science, 5(1), pp. 1–4, 1958. A. Langevin, P. Mbaraga, and J. F. Campbell. Continuous approximation models in freight distribution: An overview. Transportation Research, 30B(3):163–188, 1996. G. Laporte. Location-routing problems. In B. L. Golden and A. A. Assad, editors, Vehicle Routing: Methods and Studies, chapter 9, pages 163–197. North-Holland, Amsterdam, Netherlands, 1988.

37

G. Laporte, F. Louveaux, and H. Mercure. Models and exact solutions for a class of stochastic location-routing problems. European Journal of Operational Research, 39(1):71–78, 1989. C.-G. Lee, M. A. Epelman, C. C. White, and Y. A. Bozer. A shortest path approach to the multiple-vehicle routing problem with split pick-ups. Computers and Operations Research, 40(4):265–284, 2006. A. Mingozzi, S. Giorgi, and R. Baldacci. An exact method for the vehicle routing problem with backhauls. Transportation Science, 33(3):315–329, 1999. P. A. Mullaseril, M. Dror, and J. Leung. Split-delivery routing heuristics in livestock feed distribution. The Journal of the Operational Research Society, 48(2):107–116, 1997. G. F. Newell. Scheduling, location, transportation, and continuum mechanics; some simple approximations to optimization problems. SIAM Journal on Applied Mathematics, 25(3):346–360, 1973. I. H. Osman and N. A. Wassan. A reactive tabu search meta-heuristic for the vehicle routing problem with backhauls. Journal of Scheduling, 5(4):263–285, 2002. Y. Ouyang and C. F. Daganzo. Discretization and validation of the continuum approximation scheme for terminal system design. Transportation Science, 40(1):89–98, 2006. J.-Y. Potvin, C. Duhamel, and F. Guertin. A genetic algorithm for vehicle routing problem with backhauling. Applied Intelligence, 6:345–355, 1996. S. R. Thangiah, J.-Y. Potvin, and T. Sun. Heuristic approaches to vehicle routing with backhauls and time windows. Computers and Operations Research, 23(11):1043–1057, 1996. P. Toth and D. Vigo. An exact algorithm for the vehicle routing problem with backhauls. Transportation Science, 31:372–385, 1997. P. Toth and D. Vigo. A heuristic algorithm for the symmetric and asymmetric vehicle routing problems with backhauls. European Journal of Operational Research, 113:528–543, 1999. C. J. Vidal and M. Goetschalckx. Strategic production-distribution models: A critical review with emphasis on global supply chain models. European Journal of Operational Research, 98(1):1–18, 1997.

8

Appendix
P[λi,1,j < λi,0,j ]

For a given origin i and destination j, we show how to calculate the probability (91)

and the conditional probability P [λi,1,j > α | λi,1,j < λi,0,j ] for α ∈ (0, λi,0,j ), (92)

which is the key to completing the calculation of the expected linehaul distance E[Λi,j ] in (35). Recall that the N ? 1 terminals besides the center terminal are independent and identically distributed in a region. To facilitate the calculation of E[Λi,j ], we use the uniform distribution on a rectangular region [?a, a] × [?b, b] ? R2 , with a > 0 and b > 0. Speci?cally, terminal n = 1, 2, . . . , N ? 1 has coordinates Un ∈ R2 , where Un is uniformly distributed in [?a, a] × [?b, b], and U1 , . . . , UN ?1 are independent. The center terminal is located at (0, 0).

38

We use the L1 -norm (x, y)
1

:= |x| + |y|

and the associated metric, as in Hall (1984), Daganzo (1987), Hall (1987), Campbell (1990b), and Campbell (1993a). The L1 -norm is used for simplicity; for example, calculations using the L1 -norm are much simpler than when using the L2 -norm (x, y) 2 := x2 + y 2 . Furthermore, the ratios between the L1 and L2 -metrics between pairs of points do not vary much. To illustrate this, we show the results of the following simple experiment. We generate 10,000 pairs of points independently and uniformly distributed in the unit square, (x1 (i), y1 (i)), (x2 (i), y2 (i))), i ∈ {1, . . . , 10, 000}, and calculate the L1 and L2 -metrics between every pair of points. Figure 7 shows the cumulative distribution function of √ ratio (x2 (i), y2 (i)) ? (x1 (i), y1 (i)) 2 / (x2 (i), y2 (i)) ? the (x1 (i), y1 (i)) 1 . Note that the ratio varies between 1/ 2 and 1, and that 0.8 is a pretty good overall approximation for the ratio.
1 0.9 Cumulative fraction of points 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 Ratio of L2 to L1 norm 0.8 0.9 1

Figure 7: Cumulative distribution function of ratio ( x2 , y2 1 )?1 x1 , y1 2 .
Let (x1 , y1 ) denote the coordinates of origin i, and let (x2 , y2 ) denote the coordinates of destination j. Then, for any terminal n, λi,n,j := Un ? (x1 , y1 ) 1 + (x2 , y2 ) ? Un 1 The triangle inequality holds for any terminal n, i.e., λi,n,j ≥ (x2 , y2 ) ? (x1 , y1 )
1

= |x2 ? x1 | + |y2 ? y1 |

and λi,n,j = |x2 ? x1 | + |y2 ? y1 | when terminal n is located within the rectangle with sides parallel with the coordinate axes and having opposite corners (x1 , y1 ) and (x2 , y2 ). Note that travel through the center terminal results in a total distance λi,0,j = |x1 | + |x2 | + |y1 | + |y2 |. Thus, the random variable Λi,j de?ned in (32) satis?es |x2 ? x1 | + |y2 ? y1 | ≤ Λi,j ≤ |x1 | + |x2 | + |y1 | + |y2 |.

To calculate the probabilities in (91) and (92), we consider the following three cases separately: 1. Origin i and destination j are located in opposite quadrants. 2. Origin i and destination j are located in the same quadrant. 3. Origin i and destination j are located in adjacent quadrants. If origin i and destination j are located in opposite quadrants, then |x2 ? x1 | + |y2 ? y1 | = |x1 | + |x2 | + |y1 | + |y2 |

39

and thus Λi,j = λi,0,j = |x1 | + |x2 | + |y1 | + |y2 | for all selections of Ni and Nj . Thus E[Λi,j ] = |x1 | + |x2 | + |y1 | + |y2 |. It remains to consider the cases in which origin i and destination j are not located in opposite quadrants. Note that since the distribution of Λi,j has support on [|x2 ? x1 | + |y2 ? y1 |, |x1 | + |x2 | + |y1 | + |y2 |], it follows that the integral in the expression for E [Λi,j ] in (35) can be written as follows:
λi,0,j 0

P λi,1,j > α λi,1,j < λi,0,j

dα P λi,1,j > α λi,1,j < λi,0,j dα (93)

= |x2 ? x1 | + |y2 ? y1 | +

|x1 |+|x2 |+|y1 |+|y2 | |x2 ?x1 |+|y2 ?y1 |

Note that, because terminals 1, 2, . . . , N ? 1 are uniformly distributed in the rectangle [?a, a] × [?b, b], the probabilities in (91) and (92) can be expressed as P [λi,1,j < λi,0,j ] Area ({(x, y) ∈ [?a, a] × [?b, b] : = and P λi,1,j > α λi,1,j < λi,0,j = Area ({(x, y) ∈ [?a, a] × [?b, b] : α < (x, y) ? (x1 , y1 ) 1 + (x2 , y2 ) ? (x, y) 1 < |x1 | + |x2 | + |y1 | + |y2 |}) (95) Area ({(x, y) ∈ [?a, a] × [?b, b] : (x, y) ? (x1 , y1 ) 1 + (x2 , y2 ) ? (x, y) 1 < |x1 | + |x2 | + |y1 | + |y2 |})

(x, y) ? (x1 , y1 ) 1 + (x2 , y2 ) ? (x, y) 4ab

1

< |x1 | + |x2 | + |y1 | + |y2 |}) (94)

for α ∈ (0, λi,0,j ). The areas in the right side of (94) and (95) are straightforward to calculate. Figures 8 and 9 show the regions, with origin i and destination j located in the same quadrant and in adjacent quadrants respectively, in which location of terminal n provides a linehaul distance λi,n,j that satis?es λi,n,j ≤ α for α ∈ (|x2 ?x1 |+|y2 ?y1 |, |x1 |+|x2 |+|y1 |+|y2 |). In particular, if terminal n is located on the boundary of the eightsided region, then λi,n,j = α, and if terminal n is located in the interior of the eight-sided region, then λi,n,j < α. Setting α = λi,0,j = |x1 | + |x2 | + |y1 | + |y2 | gives the region in which location of a terminal provides a smaller linehaul distance than through the center terminal. We show that P [λi,1,j > α | λi,1,j < λi,0,j ] is a piecewise λ polynomial in α of degree at most two, from which the calculation of 0 i,0,j P [λi,1,j > α | λi,1,j < λi,0,j ] dα follows in a straightforward manner (although it is numerically dangerous due to catastrophic cancellation). The remainder of the appendix is organized as follows. We provide the detailed calculations for the case in which the rectangle [?a, a] × [?b, b] contains (x1 , y1 ) and (x2 , y2 ). The calculations for the case in which the rectangle does not contain (x1 , y1 ) and (x2 , y2 ) are similar, but involve more subcases. Section 8.1 shows the calculations for the case in which the origin and destination are located in the same quadrant, and Section 8.2 for the case in which they are located in adjacent quadrants.

40

8.1

Origin and Destination Located in the Same Quadrant

First, suppose that 0 ≤ x1 ≤ x2 and 0 ≤ y1 ≤ y2 . Then i and j are both in quadrant I. For α ∈ (|x2 ? x1 | + |y2 ? y1 |, |x1 | + |x2 | + |y1 | + |y2 |), let
S R1 (α) S R2 (α) S R3 (α) S R4 (α) S R5 (α) S R6 (α) S R7 (α) S R8 (α)

:= := := := := := := :=

Area Area Area Area Area Area Area Area

(x, y) : ?a ≤ x ≤ a and ? b ≤ y ≤ (x, y) :

?α ? x1 + x2 + y1 + y2 2
α+x1 ?x2 +y1 +y2 2

?a ≤ x ≤ ?α+x1 +x2 ?y1 +y2 2 and ?α?x1 +x2 +y1 +y2 ≤ y ≤ 2

(x, y) : ?a ≤ x ≤ a and (x, y) : (x, y) : (x, y) : (x, y) : (x, y) :

α + x1 ? x2 + y1 + y2 ≤y≤b 2 a ≤y≤
α+x1 ?x2 +y1 +y2 2 α+x1 +x2 +y1 ?y2 2

α+x1 +x2 +y1 ?y2 ≤x≤ 2 and ?α?x1 +x2 +y1 +y2 2

2(x ? y) ≥ α + x1 + x2 ? y1 ? y2 and x2 ≤ x ≤ and ?α?x1 +x2 +y1 +y2 ≤ y ≤ y1 2 ?2(x + y) ≥ α ? x1 ? x2 ? y1 ? y2 and and ?α?x1 +x2 +y1 +y2 ≤ y ≤ y1 2 2(y ? x) ≥ α ? x1 ? x2 + y1 + y2 and 2 and y2 ≤ y ≤ α+x1 ?x2 +y1 +y2

?α+x1 +x2 ?y1 +y2 2

≤ x ≤ x1

?α+x1 +x2 ?y1 +y2 2

≤ x ≤ x1

2(x + y) ≥ α + x1 + x2 + y1 + y2 and x2 ≤ x ≤ 2 and y2 ≤ y ≤ α+x1 ?x2 +y1 +y2

α+x1 +x2 +y1 ?y2 2

.

S S Figure 8 shows the regions de?ning R1 (α), . . . , R8 (α) for a value of α ∈ (|x2 ?x1 |+|y2 ?y1 |, |x1 |+|x2 |+|y1 |+|y2 |). S S The lines 1–8 in Figure 8 form the boundaries of the regions de?ning R1 (α), . . . , R8 (α), and it will be convenient to refer to the line numbers in the ?gure. For the special case with α = λi,0,j = |x1 | + |x2 | + |y1 | + |y2 |), let S Rk S := Rk (λi,0,j ).

Then, the probabilities in (94) and (95) are given by P [λi,1,j < λi,0,j ] = and P [λi,1,j > α | λi,1,j < λi,0,j ] = 4ab ?
8 k=1 S Rk

4ab
8 k=1 S Rk (α) ? 8 k=1 8 S k=1 Rk S Rk

4ab ?

(96)

Next, we allow i and j to be in any quadrant (as long as i and j are in the same quadrant), and we allow x1 > x2 or y1 > y2 . Using Figure 8, it is easy to observe that regardless of the quadrant in which i and j are, S S the areas R1 (α), . . . , R8 (α) depend only on the minimum and maximum of the distances of i and j from the x and y axes. Hence, let x := |x1 | ∧ |x2 |, x := |x1 | ∨ |x2 |, y := |y1 | ∧ |y2 |, and y := |y1 | ∨ |y2 | (97)

41

RS (α) 3

b
2 2 (x1 , α+x1 ?x2 +yr1 +y2 ) (x2 , α+x1 ?x2 +y1 +y2 ) r

? ( ?α+x1 +x2 ?y1 +y2 , y2 ) ? r 2

RS (α) ? 7 6 ? ?

?

5

d d

d

RS (α) 8
d dr(
α+x1 +x2 +y1 ?y2 , y2 ) 2 S

4d

RS (α) 2

7

(x2 , y2 ) (x1 , y1 )
d d r

r

( ?α+x1 +x2 ?y1 +y2 , y1 ) d 2 d 8 ?a

r

(x1 , ?α?x1 +x2 +y1 +y2 )(x2 , ?α?x1 +x2 +y1 +y2 ) 2 2 RS (α) 1

RS (α) d 6 dr

1

? r

? RS (α) 5

r α+x1 +x2 +y1 ?y2 , y1 ) ?( 2 ? 2? ? a

3 R4 (α)

?b Figure 8: Example of the region in which location of a terminal provides a smaller linehaul distance than α ∈ (|x2 ? x1 | + |y2 ? y1 |, |x1 | + |x2 | + |y1 | + |y2 |) if (x1 , y1 ) and (x2 , y2 ) are in the same quadrant.
Then, in general, with i and j in the same quadrant, it follows that ?α ? x + x + y + y S + b I{α<2b?x+x+y+y} R1 (α) = 2a 2 ?α + x + x ? y + y α+x?x+y+y S +a b∧ R2 (α) = 2 2 ×I{α<2a+x+x?y+y}
S R3 (α) S R4 (α)

? (?b) ∨

?α ? x + x + y + y 2

=

α+x?x+y+y 2 α+x+x+y?y = a? 2 ×I{α<2a?x?x?y+y} 2a b ? = α+x+x+y?y 1 a∧ 2 2 ×I{α<2(a+b)?x?x+y+y} 1 2 x∧ b+

I{α<2b?x+x?y?y} b∧ α+x?x+y+y 2 ? (?b) ∨ ?α ? x + x + y + y 2
2

S R5 (α)

? x∨

α+x+x?y?y ?b 2

S R6 (α) S R7 (α)

= =

?α + x + x + y + y 2

? (?a) ∨ ? (?a) ∨

?α + x + x ? y + y 2 ?α + x + x ? y + y 2
2

2

α?x?x+y+y 1 x1 ∧ b ? 2 2 ×I{α<2(a+b)+x+x?y?y} α+x+x+y?y 1 a∧ 2 2 ×I{α<2(a+b)?x?x?y?y} . ? x∨

2

S R8 (α)

=

α+x+x+y+y ?b 2

42

Note that in these area calculations one must keep in mind that part of the boundary of the eight-sided region in Figure 8 may be outside the rectangle [?a, a] × [?b, b]. Thus, depending on the boundary of the eight-sided S S S S region, some of the areas R1 (α), . . . , R8 (α) may be 0. In the expressions for R1 (α), . . . , R8 (α) above, this S is accomplished with the indicator functions (except for R6 (α), because a > 0, b > 0, and x + y > 0 imply S R6 (α) > 0). To ?nish the expected linehaul distance calculations, it follows from (93) and (96) that we must calculate
|x1 |+|x2 |+|y1 |+|y2 | |x2 ?x1 |+|y2 ?y1 |

P [λi,1,j > α | λi,1,j < λi,0,j ] dα
8 k=1 S Rk (α) ? 8 k=1 S Rk S Rk

|x1 |+|x2 |+|y1 |+|y2 |

=
|x2 ?x1 |+|y2 ?y1 | 8

4ab ?

8 k=1

dα.

(98)

S To calculate (98), we write k=1 Rk (α) as a piecewise polynomial in α of degree at most two. We will choose A(α), B(α), C1 (α), and C2 such that

A(α)α2 + B(α)α + C1 (α) =

8 S Rk (α) k=1

8

and C2 =
k=1

S Rk .

S S Considering the interval over which each indicator function in the expressions for R1 (α), . . . , R8 (α), is 0, the interval [|x2 ? x1 | + |y2 ? y1 |, |x1 | + |x2 | + |y1 | + |y2 |] can be partitioned into intervals over which A(α), B(α), and C1 (α) are constant. Also, for C(α) := C1 (α) ? C2 , 8 S Rk (α) ? k=1 k=1 8 S Rk

= =

A(α)α2 + B(α)α + C(α) ! A(α)α2 1! 2! 3!
1

[B(α)α] 2 [C(α)]

3

1, 2, 3 : 1+ 2+ 3=

=
1, 2, 3 : 1+ 2+ 3=

! A(α) 1 B(α) 2 C(α) 3 α2 ! 2! 3! 1

1+ 2

,

and thus it follows from (98) that
|x1 |+|x2 |+|y1 |+|y2 | |x2 ?x1 |+|y2 ?y1 |

P [λi,1,j > α | λi,1,j < λi,0,j ] dα
|x1 |+|x2 |+|y1 |+|y2 |

=

1 4ab ?
8 k=1 S Rk

|x2 ?x1 |+|y2 ?y1 |

1, 2, 3 : 1+ 2+ 3=

! A(α) 1 B(α) 2 C(α) 3 α2 1! 2! 3!

1+ 2

dα.

(99)

Calculating (99) is easy when the shaded region in Figure 8 is inside the rectangle [?a, a] × [?b, b] for all α ∈ (|x2 ? x1 | + |y2 ? y1 |, |x1 | + |x2 | + |y1 | + |y2 |). In this case, A(α), B(α), and C1 (α) are constant for all α ∈ (|x2 ? x1 | + |y2 ? y1 |, |x1 | + |x2 | + |y1 | + |y2 |). For example, suppose that a = b and (x1 , y1 ) = 1 1 a, a 4 4 and (x2 , y2 ) = 1 1 a, a . 2 2

Then, |x2 ? x1 | + |y2 ? y1 | = a/2 and |x1 | + |x2 | + |y1 | + |y2 | = 3a/2, and, for all α ∈ (a/2, 3a/2), 1 A(α) = ? , B(α) = 0, 2 and C1 (α) = 65 2 a . 16

43

Furthermore,
3a/2 a/2

8 k=1

S Rk = 25a2 /8. Suppose, for example, that

= 1. Then, (99) is easily calculated as
3a/2 a/2

P [λi,1,j > α | λi,1,j < λi,0,j ] dα

= = =

1 2 ? 25 a2 4a 8 1
7 2 8a 3a/2 a/2

A(α)α2 + B(α)α + C(α) dα 65 2 25 2 a ? a 16 8 dα

1 ? α2 + 2

2 a. 21

S S Note that when α = 3a/2, then R3 (α) = R4 (α) = 0. In particular, when α = 3a/2, the lines 5 and 3 are exactly on the boundary of the rectangle [?a, a] × [?a, a]. Consider a di?erent choice of (x1 , y1 ) and (x2 , y2 ), for example, 1 1 3 3 (x1 , y1 ) = a, a and (x2 , y2 ) = a, a . 4 4 4 4 S S Then R3 (α) = R4 (α) = 0 for α in a nonempty portion of the interval

(|x2 ? x1 | + |y2 ? y1 |, |x1 | + |x2 | + |y1 | + |y2 |) = (a, 2a) , which causes A(α), B(α), and C1 (α) to be piecewise constant with multiple pieces (rather than constant) over the S S interval (|x2 ? x1 | + |y2 ? y1 |, |x1 | + |x2 | + |y1 | + |y2 |). Speci?cally, R3 (α) = R4 (α) = 0 for all α ∈ [3a/2, 2a). Furthermore, A(α) = ?1 2 0 0 ?a
17 2 4 a 37 2 8 a

if if

a < α < 3a 2 3 a ≤ α < 2a 2

B(α) = C1 (α) and
8 k=1 S Rk = 21a2 /8. When 2a a

if a < α < 3 a 2 if 3 a ≤ α < 2a 2 if if a < α < 3a 2 , 3 a ≤ α < 2a 2

=

= 1, (99) becomes

P [λi,1,j > α | λi,1,j < λi,0,j ] dα
2a a

= = =

1 4a2 ? 21 a2 8 8 1 11 a2 13 a. 33
a

A(α)α2 + B(α)α + C(α) dα 17 2 21 2 a ? a 4 8
2a

3a/2

1 ? α2 + 2

dα +

3a/2

?aα +

37 2 21 2 a ? a 8 8

In general, integrating (99) requires partitioning the interval (|x2 ? x1 | + |y2 ? y1 |, |x1 | + |x2 | + |y1 | + |y2 |) into the subintervals over which A(α), B(α), and C(α) are constant. The intervals of α-values over which lines 1–8 are inside or outside the rectangle [?a, a]×[?b, b], determine the intervals of α-values over which the triangles S S S R5 (α), R7 (α), and R8 (α) have positive areas, and the intervals of α-values over which A(α), B(α), and C1 (α) are constant, from which the calculation of (99) follows directly. We ?rst consider the cases determined by which of the lines 1, 3, 5, and 7 in Figure 8 are outside the rectangle [?a, a] × [?b, b], and then, if required, consider further subcases determined by which triangles have positive area. Fortunately, dependencies imply that many of the 24 × 23 conceivable cases cannot occur. For example, examination of Figure 8 shows that if line 3 is inside the rectangle so that (α + x1 + x2 + y1 ? y2 )/2 < a, then α < 2a ? x1 ? x2 ? y1 + y2 ≤ 2a + x1 + x2 ? y1 + y2 , which implies that line 7 is also inside the rectangle. We show only the cases that can occur.

44

8.1.1

Calculating the piecewise constant functions A(α), B(α), and C(α)

For easy reference to Figure 8, and without loss of generality, assume that x1 = x, x2 = x, y1 = y, and y2 = y. S Note that in this case R6 (α) > 0 for all α ∈ (|x2 ? x1 | + |y2 ? y1 |, |x1 | + |x2 | + |y1 | + |y2 |). The following cases can occur: 1. Lines 1, 3, 5, and 7 are inside the rectangle, that is, α ∈ [|x2 ? x1 | + |y2 ? y1 |, (2a + x1 + x2 ? y1 + y2 ) ∧ (2a ? x1 ? x2 ? y1 + y2 ) ∧(2b ? x1 + x2 + y1 + y2 ) ∧ (2b ? x1 + x2 ? y1 ? y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |))
S S S which implies that R5 (α) > 0, R7 (α) > 0, and R8 (α) > 0. Then,

A(α)

= ?

1 2 1 1 2 2 (x2 ? x1 ) + (y2 ? y1 ) . 2 2

B(α) = 0 C1 (α) = 4ab +

2. Lines 1, 5, and 7 are inside the rectangle, but line 3 is outside, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2a ? x1 ? x2 ? y1 + y2 ), (2a + x1 + x2 ? y1 + y2 ) ∧ (2b ? x1 + x2 + y1 + y2 ) ∧ (2b ? x1 + x2 ? y1 ? y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |))
S S S which implies that R5 (α) > 0, R7 (α) > 0, and R8 (α) > 0. Then,

A(α)

= ?

1 4

B(α) = C1 (α)

1 (x1 + x2 ) ? a 2

1 1 1 = a(a + 4b ? x1 ? x2 ) + (x2 + x2 ) + (x2 ? x1 )2 + (y2 ? y1 )2 2 2 1 4 4

3. Lines 1, 3, and 7 are inside the rectangle, but line 5 is outside, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2b ? x1 + x2 ? y1 ? y2 ), (2a + x1 + x2 ? y1 + y2 ) ∧ (2a ? x1 ? x2 ? y1 + y2 ) ∧ (2b ? x1 + x2 + y1 + y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |))
S S S which implies that R5 (α) > 0, R7 (α) > 0, and R8 (α) > 0. Then,

A(α) B(α) C1 (α)

= ? =

1 4

1 (y1 + y2 ) ? b 2

1 1 1 2 2 = b(4a + b ? y1 ? y2 ) + (x2 ? x1 )2 + (y2 ? y1 )2 + (y1 + y2 ) 4 4 2

4. Lines 1 and 7 are inside the rectangle, but lines 3 and 5 are outside, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2a ? x1 ? x2 ? y1 + y2 ) ∨ (2b ? x1 + x2 ? y1 ? y2 ), (2a + x1 + x2 ? y1 + y2 ) ∧ (2b ? x1 + x2 + y1 + y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |))
S S which implies that R5 (α) > 0 and R7 (α) > 0. The following two subcases can occur:

(a) Part of line 4 is inside the rectangle, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2a ? x1 ? x2 ? y1 + y2 ) ∨ (2b ? x1 + x2 ? y1 ? y2 ),

(2a + x1 + x2 ? y1 + y2 ) ∧ (2b ? x1 + x2 + y1 + y2 ) ∧ (2(a + b) ? x1 ? x2 ? y1 ? y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |))

45

S which implies that R8 (α) > 0. Then,

A(α)

=

0

1 B(α) = ?(a + b) + (x1 + x2 + y1 + y2 ) 2 C1 (α) = (a + b)2 + 2ab ? a(x1 + x2 ) ? b(y1 + y2 ) + 1 2 2 2 x + x2 + y1 + y2 . 2 2 1

(b) Line 4 is outside the rectangle, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2a ? x1 ? x2 ? y1 + y2 ) ∨ (2b ? x1 + x2 ? y1 ? y2 ) ∨ (2(a + b) ? x1 ? x2 ? y1 ? y2 ), (2a + x1 + x2 ? y1 + y2 ) ∧ (2b ? x1 + x2 + y1 + y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |))
S which implies that R8 (α) = 0. Then,

1 8 1 1 B(α) = ? (a + b) + (x1 + x2 + y1 + y2 ) 2 4 1 1 2 (a + b) + 2ab + (b ? a)(x1 + x2 ? y1 ? y2 ). C1 (α) = 2 2 A(α) = ? 5. Lines 1 and 5 are inside the rectangle, but lines 3 and 7 are outside, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2a + x1 + x2 ? y1 + y2 ) ∨ (2a ? x1 ? x2 ? y1 + y2 ), (2b ? x1 + x2 + y1 + y2 ) ∧ (2b ? x1 + x2 ? y1 ? y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |))
S S S which implies that R5 (α) > 0, R7 (α) > 0, and R8 (α) > 0. Then,

A(α) C1 (α)

= =

0 2a(a + 2b) + x2 + x2 . 1 2

B(α) = ?2a

6. Lines 3 and 7 are inside the rectangle, but lines 1 and 5 are outside, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2b ? x1 + x2 + y1 + y2 ) ∨ (2b ? x1 + x2 ? y1 ? y2 ), (2a + x1 + x2 ? y1 + y2 ) ∧ (2a ? x1 ? x2 ? y1 + y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |))
S S S which implies that R5 (α) > 0, R7 (α) > 0, and R8 (α) > 0. Then,

A(α) B(α) C1 (α)

= =

0
2 2 2b(2a + b) + y1 + y2 .

= ?2b

7. Line 7 is inside the rectangle, but lines 1, 3, and 5 are outside, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2a ? x1 ? x2 ? y1 + y2 ) ∨ (2b ? x1 + x2 + y1 + y2 ) ∨ (2b ? x1 + x2 ? y1 ? y2 ), (2a + x1 + x2 ? y1 + y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |))
S which implies that R7 (α) > 0. The following three subcases can occur:

(a) Parts of lines 2 and 4 are inside the rectangle, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2a ? x1 ? x2 ? y1 + y2 ) ∨ (2b ? x1 + x2 + y1 + y2 ) ∨ (2b ? x1 + x2 ? y1 ? y2 ), (2a + x1 + x2 ? y1 + y2 ) ∧ (2(a + b) ? x1 ? x2 ? y1 ? y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |))

46

S S which implies that R5 (α) > 0 and R8 (α) > 0. Then,

A(α)

=

1 4

3 1 1 B(α) = ? a ? 2b + (x1 + x2 ) ? (y1 + y2 ) 2 2 4 1 3 2 2 C1 (α) = a2 + 4ab + 2b2 ? a(x1 + x2 ) + (x1 + x2 ) + (y1 + y2 ) ? y1 y2 . 4 4 (b) Part of line 2 is inside the rectangle and line 4 is outside the rectangle, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2a ? x1 ? x2 ? y1 + y2 ) ∨ (2b ? x1 + x2 + y1 + y2 ) ∨(2b ? x1 + x2 ? y1 ? y2 ) ∨ (2(a + b) ? x1 ? x2 ? y1 ? y2 ), (2a + x1 + x2 ? y1 + y2 ) ∧ (2(a + b) ? x1 ? x2 + y1 + y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |))
S S which implies that R5 (α) > 0 and R8 (α) = 0. Then,

A(α) B(α) C1 (α)

=

1 8

3 1 1 = ? a ? b + (x1 + x2 ? y1 ? y2 ) 2 2 4 = 3 1 2 2 (a + b)2 + 2ab ? a(x1 + x2 ) + (a + b)(y1 + y2 ) ? (x2 + x2 ) ? (y1 + y2 ) 1 2 4 4 +(x1 + x2 ? y1 ? y2 )2 .

(c) Lines 2 and 4 are outside the rectangle, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2a ? x1 ? x2 ? y1 + y2 ) ∨ (2b ? x1 + x2 + y1 + y2 ) ∨(2b ? x1 + x2 ? y1 ? y2 ) ∨ (2(a + b) ? x1 ? x2 + y1 + y2 ), (2a + x1 + x2 ? y1 + y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |))
S S which implies that R5 (α) = R8 (α) = 0. Then,

A(α)

=

0

B(α) = ?b C1 (α) 1 2 2 = b(2a + b + x1 + x2 ) + (y1 + y2 ). 2

8. Line 1 is inside the rectangle, but lines 3, 5 and 7 are outside, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2a + x1 + x2 ? y1 + y2 ) ∨ (2a ? x1 ? x2 ? y1 + y2 ) ∨ (2b ? x1 + x2 ? y1 ? y2 ), (2b ? x1 + x2 + y1 + y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |))
S which implies that R5 (α) > 0. The following three subcases can occur:

(a) Parts of lines 4 and 6 are inside the rectangle, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2a + x1 + x2 ? y1 + y2 ) ∨ (2a ? x1 ? x2 ? y1 + y2 ) ∨ (2b ? x1 + x2 ? y1 ? y2 ), (2b ? x1 + x2 + y1 + y2 ) ∧ (2(a + b) ? x1 ? x2 ? y1 ? y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |))
S S which implies that R7 (α) > 0 and R8 (α) > 0. Then,

A(α)

=

1 4

1 B(α) = ?2a ? b + (y1 + y2 ) 2 C1 (α) = 2a2 + 4ab + b2 ? b(y1 + y2 ) + 3 1 1 2 2 (x1 + x2 ) + (y1 + y2 ) + x1 x2 . 4 4 2

47

(b) Part of line 6 is inside the rectangle and line 4 is outside the rectangle, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2a + x1 + x2 ? y1 + y2 ) ∨ (2a ? x1 ? x2 ? y1 + y2 ) ∨(2b ? x1 + x2 ? y1 ? y2 ) ∨ (2(a + b) ? x1 ? x2 ? y1 ? y2 ), (2b ? x1 + x2 + y1 + y2 ) ∧ (2(a + b) + x1 + x2 ? y1 ? y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |))
S S which implies that R7 (α) > 0 and R8 (α) = 0. Then,

A(α)

=

1 8 1 (?x1 ? x2 + y1 + y2 ) 4 1 2 1 1 b + (a + b)(x1 + x2 ) + (a ? b)(y1 + y2 ) 2 2 2 1 2 + (x1 + x2 ? y1 ? y2 ) . 8

3 1 B(α) = ? a ? b + 2 2 3 2 a + 3ab + C1 (α) = 2 1 + x2 + x2 2 2 1

(c) Lines 4 and 6 are outside the rectangle, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2a + x1 + x2 ? y1 + y2 ) ∨ (2a ? x1 ? x2 ? y1 + y2 ) ∨(2b ? x1 + x2 ? y1 ? y2 ) ∨ (2(a + b) + x1 + x2 ? y1 ? y2 ), (2b ? x1 + x2 + y1 + y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |))
S S which implies that R7 (α) = R8 (α) = 0. Then,

A(α)

=

0

B(α) = ?a C1 (α) 1 = a(a + 2b + y1 + y2 ) + (x2 + x2 ). 2 2 1

9. Lines 1, 3, 5, and 7 are outside the rectangle, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2a + x1 + x2 ? y1 + y2 ) ∨ (2a ? x1 ? x2 ? y1 + y2 ) ∨(2b ? x1 + x2 + y1 + y2 ) ∨ (2b ? x1 + x2 ? y1 ? y2 ), |x1 | + |x2 | + |y1 | + |y2 |) The following ?ve subcases can occur: (a) Parts of lines 2, 4, and 6 are inside the rectangle, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2a + x1 + x2 ? y1 + y2 ) ∨ (2a ? x1 ? x2 ? y1 + y2 ) ∨(2b ? x1 + x2 + y1 + y2 ) ∨ (2b ? x1 + x2 ? y1 ? y2 ), (2(a + b) ? x1 ? x2 ? y1 ? y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |))
S S S which implies that R5 (α) > 0, R7 (α) > 0, and R8 (α) > 0. Then,

1 2 B(α) = ?2(a + b) A(α) = C1 (α) = 2(a + b)2 + 1 2 2 2 x + x2 + y1 + y2 + x1 x2 + y1 y2 . 2 2 1

(b) Parts of lines 2 and 6 are inside the rectangle and line 4 is outside the rectangle, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2a + x1 + x2 ? y1 + y2 ) ∨ (2a ? x1 ? x2 ? y1 + y2 ) ∨(2b ? x1 + x2 + y1 + y2 ) ∨ (2b ? x1 + x2 ? y1 ? y2 ) ∨ (2(a + b) ? x1 ? x2 ? y1 ? y2 ), (2(a + b) + x1 + x2 ? y1 ? y2 ) ∧ (2(a + b) ? x1 ? x2 + y1 + y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |))

48

S S S which implies that R5 (α) > 0, R7 (α) > 0, and R8 (α) = 0. Then,

A(α)

=

3 8

3 3 1 B(α) = ? a ? b ? (x1 + x2 + y1 + y2 ) 2 2 4 3 3 (a + b) + (x1 + x2 + y1 + y2 ) . C1 (α) = 2 4 (c) Part of line 6 is inside the rectangle and lines 2 and 4 are outside the rectangle, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2a + x1 + x2 ? y1 + y2 ) ∨ (2a ? x1 ? x2 ? y1 + y2 ) ∨(2b ? x1 + x2 + y1 + y2 ) ∨ (2b ? x1 + x2 ? y1 ? y2 ) ∨ (2(a + b) ? x1 ? x2 + y1 + y2 ), (2(a + b) + x1 + x2 ? y1 ? y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |))
S S S which implies that R7 (α) > 0 and R5 (α) = R8 (α) = 0. Then,

A(α)

=

1 4
2

1 1 B(α) = ?a ? b ? x1 ? x2 2 2 1 C1 (α) = a + b + (x1 + x2 ) 2

+

1 2 (y1 + y2 ) . 4

(d) Part of line 2 is inside the rectangle and lines 4 and 6 are outside the rectangle, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2a + x1 + x2 ? y1 + y2 ) ∨ (2a ? x1 ? x2 ? y1 + y2 ) ∨(2b ? x1 + x2 + y1 + y2 ) ∨ (2b ? x1 + x2 ? y1 ? y2 ) ∨ (2(a + b) + x1 + x2 ? y1 ? y2 ), (2(a + b) ? x1 ? x2 + y1 + y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |))
S S S which implies that R5 (α) > 0 and R7 (α) = R8 (α) = 0. Then,

A(α)

=

1 4
2

1 1 B(α) = ?a ? b ? y1 ? y2 2 2 1 C1 (α) = a + b + (y1 + y2 ) 2 (e) Lines 2, 4, and 6 are outside the rectangle, that is, α

+

1 2 (x1 + x2 ) . 4

∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2a + x1 + x2 ? y1 + y2 ) ∨ (2a ? x1 ? x2 ? y1 + y2 ) ∨(2b ? x1 + x2 + y1 + y2 ) ∨ (2b ? x1 + x2 ? y1 ? y2 ) ∨ (2(a + b) + x1 + x2 ? y1 ? y2 ) ∨(2(a + b) ? x1 ? x2 + y1 + y2 ), |x1 | + |x2 | + |y1 | + |y2 |)

S S S which implies that R5 (α) = R7 (α) = R8 (α) = 0. Then,

A(α)

=

1 8

1 1 B(α) = ? (a + b) ? (x1 + x2 + y1 + y2 ) 2 4 1 (a + b) (a + b + x1 + x2 + y1 + y2 ) C1 (α) = 2 1 2 + (x1 + x2 + y1 + y2 ) . 8

49

8.2

First, suppose that x1 ≤ 0 ≤ x2 and 0 ≤ y1 ≤ y2 . Then i is in quadrant II and j is in quadrant I, as in Figure 9. For α ∈ (|x2 ? x1 | + |y2 ? y1 |, |x1 | + |x2 | + |y1 | + |y2 |), let
A R1 (α) A R2 (α) A R3 (α) A R4 (α) A R5 (α) A R6 (α) A R7 (α) A R8 (α)

:= := := := := := := :=

Area Area Area Area Area Area Area Area

(x, y) : ?a ≤ x ≤ a and ? b ≤ y ≤ (x, y) :

?α ? x1 + x2 + y1 + y2 2
α+x1 ?x2 +y1 +y2 2

?a ≤ x ≤ ?α+x1 +x2 ?y1 +y2 2 and ?α?x1 +x2 +y1 +y2 ≤ y ≤ 2

(x, y) : ?a ≤ x ≤ a and (x, y) : (x, y) : (x, y) : (x, y) : (x, y) :

α + x1 ? x2 + y1 + y2 ≤y≤b 2 a ≤y≤
α+x1 ?x2 +y1 +y2 2 α+x1 +x2 +y1 ?y2 2

α+x1 +x2 +y1 ?y2 ≤x≤ 2 and ?α?x1 +x2 +y1 +y2 2

2(x ? y) ≥ α + x1 + x2 ? y1 ? y2 and x2 ≤ x ≤ and ?α?x1 +x2 +y1 +y2 ≤ y ≤ y1 2 2(x + y) ≥ ?α + x1 + x2 + y1 + y2 and and ?α?x1 +x2 +y1 +y2 ≤ y ≤ y1 2 2(y ? x) ≥ α ? x1 ? x2 + y1 + y2 and 2 and y2 ≤ y ≤ α+x1 ?x2 +y1 +y2

?α+x1 +x2 ?y1 +y2 2

≤ x ≤ x1

?α+x1 +x2 ?y1 +y2 2

≤ x ≤ x1

2(x + y) ≥ α + x1 + x2 + y1 + y2 and x2 ≤ x ≤ 2 and y2 ≤ y ≤ α+x1 ?x2 +y1 +y2

α+x1 +x2 +y1 ?y2 2

.

A A Figure 9 shows the regions de?ning R1 (α), . . . , R8 (α) for a value of α ∈ (|x2 ?x1 |+|y2 ?y1 |, |x1 |+|x2 |+|y1 |+|y2 |). A A The lines 1–8 in Figure 9 form the boundaries of the regions de?ning R1 (α), . . . , R8 (α). For the special case with α = λi,0,j = |x1 | + |x2 | + |y1 | + |y2 |), let A A Rk := Rk (λi,0,j ).

Then, the probabilities in (94) and (95) are given by P [λi,1,j < λi,0,j ] = and P [λi,1,j > α | λi,1,j < λi,0,j ] = 4ab ?
8 k=1 A Rk

4ab
8 k=1 A Rk (α) ? 8 k=1 8 A Rk k=1 A Rk

4ab ?

(100)

Next, we allow i and j to be in any adjacent quadrants, and we allow x1 > x2 or y1 > y2 . It is easy to see how to choose a new coordinate system with the same origin and scale as the original coordinate system and with axes coinciding with the axes of the original coordinate system, so that j will be in the new quadrant I and i will be in the new quadrant II. Furthermore, using Figure 9, it is easy to observe that in the new coordinate A A system, the areas R1 (α), . . . , R8 (α) depend only on the minimum and maximum of the distances of i and j from the new x axis. Hence, if x1 and x2 have opposite signs and y1 and y2 have the same sign, then let x := ?|x1 |, x := |x2 |, y := |y1 | ∧ |y2 |, y := |y1 | ∨ |y2 |, a := a, and b := b (101)

Otherwise, if x1 and x2 have the same sign and y1 and y2 have opposite signs, then let x := ?|y1 |, x := |y2 |, y := |x1 | ∧ |x2 |, y := |x1 | ∨ |x2 |, a := b, and b := a (102)

50

b
2 (x1 , α+x1 ?x2 +y1 +y2 ) r ? RA (α) 7

RA (α) 3 4

2 (x2 , α+x1 ?x2 +y1 +y2 )

r ( x1 +x2 ?y1 +y2 ?α , y2 )? 2 6 RA (α) 2

? 5

(x2 , y2 )

r d A (α) R8 d r 3 dr α+x1 +x2 +y1 ?y2 , y2 ) (
2

2

RA (α) 4

( x1 +x2 ?y1 +y2 ?α , y1 )r (x ,ry ) 2 d 1 1 d7 RA (α) 6 dr ?a

8

(x1 , ?x1 +x2 +y1 +y2 ?α ) 2

(x2 , ?x1 +x2 +y1 +y2 ?α ) 2

r α+x1 +x1 +y1 ?y2 , y1 ) ( 2 1? ? (α) A R ? r 5

a

RA (α) 1

?b Figure 9: Example of the region in which location of a terminal provides a smaller L1 linehaul distance than α ∈ (|x2 ? x1 | + |y2 ? y1 |, |x1 | + |x2 | + |y1 | + |y2 |) if (x1 , y1 ) and (x2 , y2 ) are in adjacent quadrants.
Then
A R1 (α) A R2 (α)

= ?a α + a (?x + x + y + y) + 2a b α+x?x+y+y ?α ? x + x + y + y = b ∧ ? 2 2 ×I{α<2a +x+x?y+y} = α+x?x+y+y 2 α+x+x+y?y = a ? 2 ×I{α<2a ?x?x?y+y} 2a b ? = = = 1 2 1 2 a ∧ α+x+x+y?y 2 I{α<2b ?x+x?y?y} b ∧ α+x?x+y+y 2
2

a +

?α + x + x ? y + y 2

A R3 (α) A R4 (α)

?

?α ? x + x + y + y 2

A R5 (α) A R6 (α) A R7 (α)

?x
2

x ? (?a ) ∨

?α + x + x ? y + y 2

α?x?x+y+y 1 x∧ b ? 2 2 ×I{α<2(a +b )+x+x?y?y} α+x+x+y?y 1 a ∧ 2 2 ×I{α<2(a +b )?x?x?y?y}

? (?a ) ∨

?α + x + x ? y + y 2
2

2

A R8 (α)

=

? x∨

α+x+x+y+y ?b 2

51

A Similar to Section 8.1, we write k=1 Rk (α) as a piecewise polynomial in α of degree at most two, and we choose A(α), B(α), C1 (α), and C2 such that

8

A(α)α2 + B(α)α + C1 (α) =

8 A Rk (α) k=1

8

and C2 =
k=1

A Rk .

8.2.1

Calculating the piecewise constant functions A(α), B(α), and C(α)

For easy reference to Figure 9, and without loss of generality, assume that x1 = x, x2 = x, y1 = y, y2 = y, a = a , A A and b = b . Note that in this case R5 (α) > 0 and R6 (α) > 0 for all α ∈ (|x2 ? x1 | + |y2 ? y1 |, |x1 | + |x2 | + |y1 | + |y2 |). As in Section 8.1, when partitioning the interval (|x2 ? x1 | + |y2 ? y1 |, |x1 | + |x2 | + |y1 | + |y2 |) into the subintervals over which A(α), B(α), and C(α) are constant, we show only the cases that can occur, namely the following cases: 1. Lines 2, 4, and 6 are inside the rectangle, that is, α ∈ [|x2 ? x1 | + |y2 ? y1 |, (2a + x1 + x2 ? y1 + y2 ) ∧ (2b ? x1 + x2 ? y1 ? y2 ) ∧ (2a ? x1 ? x2 ? y1 + y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |))
A A which implies that R7 (α) > 0 and R8 (α) > 0. Then,

A(α)

= ?

1 2 1 2 2 2 x + x2 + y1 + y2 ? x1 x2 ? y1 y2 . 2 2 1

B(α) = 0 C1 (α) = 4ab +

2. Lines 2 and 4 are inside the rectangle, but line 6 is outside, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2a + x1 + x2 ? y1 + y2 ), (2b ? x1 + x2 ? y1 ? y2 ) ∧ (2a ? x1 ? x2 ? y1 + y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |))
A A which implies that R7 (α) > 0 and R8 (α) > 0. Then,

A(α)

= ?

1 4

1 B(α) = ?a ? (x1 + x2 ) 2 C1 (α) = a(a + 4b + x1 + x2 ) + 1 3 2 1 2 x + x2 + (y2 ? y1 ) ? x1 x2 . 2 4 1 4 2

3. Lines 4 and 6 are inside the rectangle, but line 2 is outside, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2a ? x1 ? x2 ? y1 + y2 ), (2a + x1 + x2 ? y1 + y2 ) ∧ (2b ? x1 + x2 ? y1 ? y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |))
A A which implies that R7 (α) > 0 and R8 (α) > 0. Then,

A(α)

= ?

1 4

1 B(α) = ?a + (x1 + x2 ) 2 C1 (α) = a(a + 4b ? x1 ? x2 ) + 1 3 2 1 2 x + x2 + (y2 ? y1 ) ? x1 x2 . 2 4 1 4 2

52

4. Lines 2 and 6 are inside the rectangle, but line 4 is outside, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2b ? x1 + x2 ? y1 ? y2 ), (2a + x1 + x2 ? y1 + y2 ) ∧ (2a ? x1 ? x2 ? y1 + y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |))
A A which implies that R7 (α) > 0 and R8 (α) > 0. Then,

A(α) B(α) C1 (α)

= ?

1 4

1 = ?b + (y1 + y2 ) 2 = b(4a + b ? y1 ? y2 ) + 1 1 3 2 2 2 (x2 ? x1 ) + y + y2 ? y1 y2 . 4 4 1 2

5. Line 2 is inside the rectangle, but lines 4 and 6 are outside, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2a + x1 + x2 ? y1 + y2 ) ∨ (2b ? x1 + x2 ? y1 ? y2 ), (2a ? x1 ? x2 ? y1 + y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |)) which implies that R8 (α) > 0. The following two subcases can occur: (a) Part of line 5 is inside the rectangle, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2a + x1 + x2 ? y1 + y2 ) ∨ (2b ? x1 + x2 ? y1 ? y2 ), (2a ? x1 ? x2 ? y1 + y2 ) ∧ (2(a + b) + x1 + x2 ? y1 ? y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |)) which implies that R7 (α) > 0. Then, A(α) = 0 1 1 B(α) = ?(a + b) ? (x1 + x2 ) + (y1 + y2 ) 2 2 C1 (α) = a2 + b2 + 4ab + a(x1 + x2 ) ? b(y1 + y2 ) + 1 2 2 2 x + x2 + y1 + y2 . 2 2 1

(b) Line 5 is outside the rectangle, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2a + x1 + x2 ? y1 + y2 ) ∨(2b ? x1 + x2 ? y1 ? y2 ) ∨ (2(a + b) + x1 + x2 ? y1 ? y2 ), (2a ? x1 ? x2 ? y1 + y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |)) which implies that R7 (α) = 0. Then, 1 8 1 1 1 B(α) = ? (a + b) ? (x1 + x2 ) + (y1 + y2 ) 2 4 4 1 1 2 2 a + b + 3ab + (a ? b)(x1 + x2 + y1 + y2 ) C1 (α) = 2 2 1 3 2 2 2 2 + x1 + x2 + y1 + y2 + (?x1 x2 + x1 y1 + x1 y2 + x2 y1 + x2 y2 ? y1 y2 ) . 8 4 A(α) = ? 6. Line 6 is inside the rectangle, but lines 2 and 4 are outside, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2a ? x1 ? x2 ? y1 + y2 ) ∨ (2b ? x1 + x2 ? y1 ? y2 ), (2a + x1 + x2 ? y1 + y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |)) which implies that R7 (α) > 0. The following two subcases can occur:

53

(a) Part of line 3 is inside the rectangle, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2a ? x1 ? x2 ? y1 + y2 ) ∨ (2b ? x1 + x2 ? y1 ? y2 ), (2a + x1 + x2 ? y1 + y2 ) ∧ (2(a + b) ? x1 ? x2 ? y1 ? y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |)) which implies that R8 (α) > 0. Then, A(α) B(α) C1 (α) = 0

1 = ?(a + b) + (x1 + x2 + y1 + y2 ) 2 = 1 2 2 4ab + a2 + b2 ? a(x1 + x2 ) ? b(y1 + y2 ) + (x2 + x2 + y1 + y2 ). 2 2 1

(b) Line 3 is outside the rectangle, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2a ? x1 ? x2 ? y1 + y2 ) ∨(2b ? x1 + x2 ? y1 ? y2 ) ∨ (2(a + b) ? x1 ? x2 ? y1 ? y2 ), (2a + x1 + x2 ? y1 + y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |)) which implies that R8 (α) = 0. Then, 1 8 1 1 1 B(α) = ? (a + b) + (x1 + x2 ) + (y1 + y2 ) 2 4 4 1 2 1 2 C1 (α) = 3ab + (a + b ) + (b ? a)(x1 + x2 ? y1 ? y2 ) 2 2 1 2 1 2 2 2 + (x1 + x2 + y1 + y2 ) ? (x1 + x2 + y1 + y2 )2 . 2 8 A(α) = ? 7. Line 4 is inside the rectangle, but lines 2 and 6 are outside, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2a ? x1 ? x2 ? y1 + y2 ) ∨ (2a + x1 + x2 ? y1 + y2 ), (2b ? x1 + x2 ? y1 ? y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |)) which implies that R7 (α) > 0 and R8 (α) > 0. Then, A(α) C1 (α) = = 0 4ab + 2a2 + x2 + x2 . 1 2

B(α) = ?2a

8. Lines 2, 4, and 6 are outside the rectangle, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2a ? x1 ? x2 ? y1 + y2 ) ∨ (2a + x1 + x2 ? y1 + y2 ) ∨ (2b ? x1 + x2 ? y1 ? y2 ), |x1 | + |x2 | + |y1 | + |y2 |) The following four subcases can occur: (a) Parts of lines 3 and 5 are inside the rectangle, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2a ? x1 ? x2 ? y1 + y2 ) ∨ (2a + x1 + x2 ? y1 + y2 ) ∨ (2b ? x1 + x2 ? y1 ? y2 ), (2(a + b) + x1 + x2 ? y1 ? y2 ) ∧ (2(a + b) ? x1 ? x2 ? y1 ? y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |))

54

which implies that R7 (α) > 0 and R8 (α) > 0. Then, A(α) = 1 4 1 (y1 + y2 ) 2 1 1 2 1 2 2 x + x2 + (x1 + x2 ) + (y1 + y2 ) 2 2 1 4 4

B(α) = ?2a ? b + C1 (α) =

2a2 + b2 + 4ab ? b (y1 + y2 ) +

(b) Part of line 5 is inside the rectangle and line 3 is outside the rectangle, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2a ? x1 ? x2 ? y1 + y2 ) ∨ (2a + x1 + x2 ? y1 + y2 ) ∨(2b ? x1 + x2 ? y1 ? y2 ) ∨ (2(a + b) ? x1 ? x2 ? y1 ? y2 ), (2(a + b) + x1 + x2 ? y1 ? y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |)) which implies that R7 (α) > 0 and R8 (α) = 0. Then, A(α) = 1 8

3 1 1 B(α) = ? a ? b ? (x1 + x2 ? y1 ? y2 ) 2 2 4 3 2 1 2 1 1 a + b + 3ab + a (x1 + x2 + y1 + y2 ) + b (x1 + x2 ? y1 ? y2 ) C1 (α) = 2 2 2 2 1 1 2 + x2 + x2 + (x1 + x2 ? y1 ? y2 ) . 2 2 1 8 (c) Part of line 3 is inside the rectangle and line 5 is outside the rectangle, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2a ? x1 ? x2 ? y1 + y2 ) ∨ (2a + x1 + x2 ? y1 + y2 ) ∨(2b ? x1 + x2 ? y1 ? y2 ) ∨ (2(a + b) + x1 + x2 ? y1 ? y2 ), (2(a + b) ? x1 ? x2 ? y1 ? y2 ) ∧ (|x1 | + |x2 | + |y1 | + |y2 |)) which implies that R7 (α) = 0 and R8 (α) > 0. Then, A(α) = 1 8

3 1 1 B(α) = ? a ? b + (x1 + x2 + y1 + y2 ) 2 2 4 3 2 1 2 1 1 a + b + 3ab ? a (x1 + x2 ? y1 ? y2 ) ? b (x1 + x2 + y1 + y2 ) C1 (α) = 2 2 2 2 1 1 2 + x2 + x2 + (x1 + x2 + y1 + y2 ) . 2 2 1 8 (d) Lines 3 and 5 are outside the rectangle, that is, α ∈ [(|x2 ? x1 | + |y2 ? y1 |) ∨ (2a ? x1 ? x2 ? y1 + y2 ) ∨ (2a + x1 + x2 ? y1 + y2 ) ∨(2b ? x1 + x2 ? y1 ? y2 ) ∨ (2(a + b) + x1 + x2 ? y1 ? y2 ) ∨ (2(a + b) ? x1 ? x2 ? y1 ? y2 ), |x1 | + |x2 | + |y1 | + |y2 |) which implies that R7 (α) = 0 and R8 (α) = 0. Then, A(α) B(α) C1 (α) = 0 1 2 x + x2 . 2 2 1 = ?a = a2 + 2ab + a(y1 + y2 ) +

55

9
9.1

Online Appendix
Computational Results for MILP formulation

For very small problem instances, the mixed integer linear program (3) can be solved with available software. Table 3 shows CPU times for some small instances obtained with CPlex 9.0. All instances were randomly generated, with a single scenario each.

Number of Origins Number of Destinations Number of Candidate Terminals Instance 1 Instance 2 Instance 3 Instance 4 Instance 5 Mean Standard Deviation

4 4 3 62.62 308.35 31.65 184.23 490.47 215.46 188.67

5 5 3 182.23 36895.92 692.00 5002.46 1873.43 8929.21 15745.65

6 6 3 43164.04 22475.01 Stop after 187643.12 Stop after 115830.61 Stop after 379172.34 N/A N/A

Table 3: CPU time in seconds for MILP formulation

9.2

Data Sets

The data sets used for the results reported in Sections 7 and 9.1 can be found at http://www.scl.gatech.edu/research/casestudies/. The distances dij between pairs of points i and j were given by the least great-circle distances in miles between the pairs of points. Also, for the purpose of the continuous approximation, (latitude, longitude) coordinates were converted to Cartesian (x, y) coordinates according to the Albers Equal Area Conic Projection Method. The L1 metric between pairs of (x, y) coordinates were multiplied with a factor of 3150, which gives approximately the least great-circle distance in miles between the pair of points. All scenarios were assigned equal weights p(ω). In addition, the following parameters were used: Fixed cost per terminal per time period cm = \$10000. Transportation cost per vehicle-mile = \$1. Cost per time period for each vehicle based at each terminal Cv = \$100. Cost for each vehicle that is used during a time period cv = \$100. Vehicle capacity used in detailed vehicle routing calculations Qv = 3000ft3 . Vehicle capacity used in continuous approximation calculations Qv = 2900ft3 . β-coe?cient for approximating detour distance β = 2. The rectangle [?a, a] × [?b, b] containing the terminals was chosen as follows: Let (xi , yi ) denote the coordinates of origin or destination i. Let x := maxi∈O∪D xi , x := mini∈O∪D xi , y := maxi∈O∪D yi , y := mini∈O∪D yi . Then a := (x ? x)/2 and b := (y ? y)/2.

56