An Improved Algorithm for Fixed-Hub Single Allocation Problem
This paper discusses the fixed-hub single allocation problem (FHSAP). In this problem, a network consists of hub nodes and terminal nodes. Hubs are fixed and fully connected; each terminal node is connected to a single hub which routes all its traffic. The goal is to minimize the cost of routing the traffic in the network. In this paper, we propose a linear programming (LP)-based rounding algorithm. The algorithm is based on two ideas. First, we modify the LP relaxation formulation introduced in Ernst and Krishnamoorthy (1996, 1999) by incorporating a set of validity constraints. Then, after obtaining a fractional solution to the LP relaxation, we make use of a geometric rounding algorithm to obtain an integral solution. We show that by incorporating the validity constraints, the strengthened LP often provides much tighter upper bounds than the previous methods with a little more computational effort, and the solution obtained often has a much smaller gap with the optimal solution. We also formulate a robust version of the FHSAP and show that it can guard against data uncertainty with little cost.
Key words: hub location; network design; linear programming; worst-case analysis
Hub-and-spoke networks have been widely used in transportation, logistics, and telecommunication systems. In such networks, traffic is routed from numerous nodes of origin to specific destinations through hub facilities. The use of hub facilities allows for the replacement of direct connections between all nodes with fewer, indirect connections. One main benefit is the economies of scale as a result of the consolidation of flows on relatively few arcs connecting the nodes. In the United States, hub-and-spoke routing is practically universal. Airlines adopted it after the industry was deregulated in . Many logistics service providers such as UPS and Fedex also have distribution systems using the hub-and-spoke structure.
Given its widespread use, it is of practical importance to design efficient hub-and-spoke networks. In the literature, such problems are often referred as the hub location problems, in which two major questions are studied: 1) where the hubs should be located and 2) how the traffic/flow should be routed. We refer the readers to  for a comprehensive review of the literature on the hub location problems.
In this paper, we focus on a sub-problem which is called the fixed-hub single allocation problem (FHSAP). In the FHSAP, the locations of the hubs are fixed, and the decisions are to assign each terminal node to a unique hub. Although the FHSAP is a sub-problem of the hub location problems, it is still of great interest. First, in many practical situations, the locations of the hubs are pre-determined and remain unchanged in a medium to long term. In such cases, the hubs can be viewed as fixed and only the assignment of the terminal nodes needs to be decided. Second, the number of nodes that can be used as hubs are usually small, which makes it possible to enumerate all possible locations of the hubs to find the optimal location. Therefore, solving the fixed-hub allocation problem efficiently would be of great help for solving the hub location problem. Moreover, even confined to fixed hubs, optimally assigning terminal nodes to hubs is still a challenging task. Indeed, it is known that FHSAP is NP-Hard even for problem with 3 hubs . Therefore, designing efficient algorithms to solve FHSAP is of great interest, both to researchers and practitioners.
To address the FHSAP, several prior approaches have been proposed. In O’Kelly , the author proposes a quadratic integer program to model this problem. The formulation is non-convex and thus hard to solve. Therefore, the author proposes two heuristic methods to solve it. Following , several other heuristic methods are proposed, see, e.g., Klincewicz , Campbell  and Skorin-Kapov et al. .
One major method to solve the FHSAP is to use a linearization model for the quadratic integer program in . Such linearizations are developed in [2, 13, 12, 6, 7]. One of the earliest such linearization model is introduced by , in which a natural LP relaxation of the quadratic integer program is obtained. This LP relaxation is quite attractive: Skorin-Kapor et al  show that it is very tight and outputs integral solutions automatically in of the instances they test. However, the size of this LP relaxation is relatively large and thus restricts its applications to large-scale problems. To solve this problem, Ernst and Krishnamoorthy [6, 7] propose a further relaxation of the model. The idea of the further relaxation is to use combined flow variables, and the size of the further relaxed LP is significantly smaller than that in  and . However, in some situations, the further relaxed model has a large gap with the optimal solution.
In this paper, we propose a new LP relaxation for the FHSAP. Our new LP relaxation is based on the one proposed in [6, 7], but we add a set of flow validity constraints to it. We show that by adding the flow validity constraints, we can often tighten the gap between the LP relaxation of [6, 7] and the optimal solution, and yield integral solutions more frequently. Moreover, it comes with reasonable computational cost. Therefore, we believe our approach is a good balance between the LP relaxation by  and [6, 7].
Besides finding a suitable LP relaxation, another important question is how to round a fractional solution of the LP into a feasible solution to the FHSAP. In this paper, we adopt a geometric rounding algorithm introduced by Ge et al. . In , the authors propose a random geometric rounding scheme for a class of assignment problems. They prove that this rounding technique can be applied to the FHSAP and lead to a constant ratio approximation algorithm under the equilateral structure. In this paper, we show that our newly proposed LP relaxation combined with the geometric rounding algorithm yields good solutions to the FHSAP efficiently. It is worth noting that another dependent rounding scheme by Kleinberg and Tardos  can also be adopted to round the solutions.
In practical cases, the demands in the FHSAP may be unknown. To tackle such situations, we propose a robust programming approach for the FHSAP when the demands are only known to be within a certain convex set. We derive a convex programming relaxation for the robust formulation which can be solved efficiently. We show in our numerical tests that by employing the decisions of the robust model, we can guard against the demand uncertainty with little cost, therefore it might be of practical interest.
The remainder of the paper is organized as follows. In Section 2, we introduce the model and the LP relaxation we propose for the FHSAP. Then we introduce the geometric rounding scheme in Section 3. In Section 4, we perform numerical tests to show that our proposed approach can indeed obtain better solutions to this problem. In Section 5, we establish a robust model for the FHSAP, and study the solution of the robust model. Then we conclude our paper in Section 6.
2 Model and Formulation
This section defines the fixed-hub single allocation problem, reviews and modifies previously proposed mathematical programs. By the terminology of communication networks, the problem is to build a two-level network consisting of hubs and terminal nodes , see Figure 1 for an illustration. In the FHSAP, we assume that there are fixed hubs denoted by (airports, routers, concentrators, etc.), which are transit nodes that are used to route traffic. There are terminals nodes denoted by (cities, computers, etc.) which represent the origins and the destinations of the traffic. Here all hubs are fully connected and each terminal node is connected to exactly one hub.
In this network, there is a demand to be routed from to , for each pair of terminal nodes and . In order to route the demands between two terminal nodes, the origin node has to deliver all its demands to the hub it is assigned to. Then this hub sends them to the hub the destination node is assigned to (this step is skipped if both nodes are assigned to the same hub). Finally the destination node gets the demands from its hub. No direct routing between two terminal nodes is permitted. Two types of costs are counted during the transportation, a per unit transportation cost to transport demand from terminal node to hub and a per unit transportation cost to transport demand from hub to hub . The problem is to assign a hub for each terminal node such that the total transportation cost is minimized.
The first mathematical formulation to study the FHSAP is by O’Kelly , in which he formulates the problem as a quadratic integer program.111In fact, his formulation is for a more general problem, the uncapacitated single allocation -hub median problem. In this paper, we only confine our discussion to the FHSAP and thus adapt his formulation (and later formulations) to the FHSAP. Define to be the assignment variables. The quadratic formulation for the FHSAP is:
Here we assume that all coefficients are non-negative and , , for all and . Note that the transportation cost from cities to hubs, , is linear in . Later we call it the linear cost and denote it by . Similarly, we call the other part of the objective function the inter-hub cost or quadratic cost, and denote it by .
Campbell  linearized O’Kelly’s model by formulating a mixed integer linear program (MILP) as follows:
Here is the portion of the flow from city to city via hub and sequentially. The formulation involves nonnegative variables and constraints. This formulation enables us to obtain an LP relaxation for the FHSAP by replacing the zero-one constraints with non-negative constraints. In the following, we refer this LP relaxation as LP1. As shown in , LP1 is usually very tight and often produces integral solutions. However, the size of LP1 is relative large, which restricts its applications to large-scale problems.
In order to reduce the solution complexity, Ernst and Krishnamoorthy [6, 7] propose a flow formulation to obtain a further relaxation of this problem. In this formulation, one does not need to specify the route for a pair of terminal nodes and , i.e., one does not need the decision variable . Instead, one defines where is the total amount of the flow originated from city and routed from hub to a different hub . Then the FHSAP can be further relaxed to:
where and denote the total demands from and to respectively. Note that this modified formulation involves only nonnegative variables and linear constraints, which decreases from that of MILP1 by a factor of . We can then obtain an LP relaxation from MILP2, which we denote by LP2.
To see that MILP2 is indeed a further relaxation of the problem, note that any feasible assignment to the FHSAP with the flow vector is always a feasible solution to MILP2 with the objective value equal to the transportation cost. Since MILP2 reduces the formulation size by , its LP relaxation is also easier to solve. However, despite that it is proved that MILP2 is an exact formulation when all the costs in the system are equal, in general, there might be a positive gap between the optimal value of MILP2 and the true optimal solution. And in our numerical tests, we find that the gap sometimes is quite large. Therefore, it is useful to find an improved formulation of MILP2 without adding too much complexity.
In the following, we propose a stronger formulation than MILP2. The main idea is to add a set of validity constraints based on the following observation.
Let and be defined as in MILP2. For any and , we have
We verify equation (2) in two cases.
If , then
If , then
Therefore, equation (2) holds in both cases. ∎
We call this problem MILP2’ and its LP relaxation LP2’. Note that LP2’ has both variables and constraints. We further reduce the number of additional constraints by summing up the validity constraints. We get our final formulation as follows:
We call the LP relaxation of MILP3 by LP3. The number of variables and constraints in MILP3 and LP3 are both . Although it doesn’t reduce the size of LP2’ significantly, computational results indicate that LP3 can be solved much more efficiently, yet the results are usually quite good.
The above formulations can serve two purposes. First, it provides a tighter lower bound for the FHSAP than LP2. Second, it provides a new way to solve the FHSAP using LP relaxations. In the next section, we show how to obtain an integral solution from the fractional solution solved from the LP relaxations. In Section 4, we perform numerical tests to show the performance of our proposed approach.
3 Rounding Procedure: A Geometric Rounding Algorithm
Note that in the above formulations, a solution to the FHSAP can be completely defined by the assignment variables . Therefore, after solving an LP relaxation (LP1, LP2 or LP3), we only need to focus on rounding the fractional assignment variables to binary integers. Note that in the three relaxations presented above (LP1, LP2 or LP3), we all have the constraints that . Therefore, for a terminal node , any optimal solution of the LP relaxation must fall on the standard dimensional simplex:
A fractional assignment vector on node corresponds to a non-vertex point of . Our goal is to round any fractional solution to a vertex point of , which is of the form:
It is clear that has exactly vertices. In the following, we denote the vertices of by , where the th coordinate of is .
Before presenting the rounding procedure, we define some notations of the geometry of the problem. For a point , connect with all vertices of . Denote the polyhedron which exactly has vertices by . Thus the simplex can be partitioned into polyhedrons , and the interiors of any distinct pair of these polyhedrons do not intersect.
We are now ready to present our randomized rounding algorithm. Note that this rounding procedure is applicable to problems other than FHSAP, as long as the feasible set of the problems is the set of vertices of a simplex.
Geometric Rounding Algorithm (FHSAP-GRA):
Solve an LP relaxation of the FHSAP (LP1, LP2 or LP3). Denote the optimal solution by .
Generate a random vector , which follows a uniform distribution on .
For each , if falls into , let ; other components .
An illustration of this procedure is shown in Figure 2. Next we discuss some theoretical properties of this rounding technique. We first consider the uniform cost case, in which the inter-hub costs are all equal. In this case, an important observation is that for any integral solution
The above observation is made by Kleinberg and Tardos  and some other literature thereafter. Furthermore, by exploring the structure of the geometric rounding, Ge et al.  proved the following two properties of this rounding technique.
For any given , .
For any and , define . Then
For any , if we randomly round and to vertices and in by the procedure in FHSAP-GRA, then
Furthermore, by Chekuri et al , if is an optimal solution to LP1, then
Assume all ’s are equal. Let be the optimal solution of LP1. And be the integral solution obtained from FHSAP-GRA. Then is a 2-approximation for the original problem. That is
where OPT is the optimal value of FHSAP-QP.
In a more general case when the inter-hub costs are not equal, we have the following extension of Theorem 3.
Let , and . Then the algorithm FHSAP-GRA using the LP relaxation LP1 has a performance guarantee of .
Given an instance of the FHSAP, we build another instance denoted by , in which all of the inter-hub costs are set to . Obviously, for any allocation vector , the cost in is less than times the cost in . And by using the algorithm FHSAP-GRA, we can obtain a 2-approximation for , which directly translates into a -approximation for . ∎
It is tempting to extend the above results to the case when LP3 is used. However, since LP3 is a further relaxation of LP1 and the structure of the objective is not the same, we are not able to prove the same bound when we use GRA combined with LP3. However, as we will show in our numerical experiments, using LP3 combined with the geometric rounding procedure produces good solutions.
4 Computational Results
In this section, we implement our algorithm (FHSAP-GRA) and report its performance. We test the algorithm on both randomly generated instances (Table 1) and a benchmark problem (Table 2). All linear programs in the experiments are solved by CPLEX version at a workstation with 3GHz CPU and 8GB memory.
In Table 1, we consider three setups of the problem: , ; , and , . In each of the setup, demands between any two cities are generated from uniform distributions and hub to city costs are generated from . Then we choose different distributions for the inter-hub costs to conduct our tests, which are shown in the second column. We try to solve all the three LP relaxation problems (LP1, LP2 and LP3) introduced in Section 2 and apply the geometric rounding algorithm to the solutions we obtain. The results are shown in the three columns GRA-LP1, GRA-LP2 and GRA-LP3 respectively. Within each of the three sets of experiments, the CPU columns show the time (in seconds) our program takes to solve each LP relaxation (we find that the time to perform the rounding procedure is negligible, therefore we only report the time to solve the LPs in our test results). Gap-LP1 columns show the gap between the cost of our obtained integral solution (for each fractional solution we obtain, we do the random rounding 5000 times and pick the best results) and the optimal value of LP1. More precisely, if we denote the optimal value of LPi by and the value of an integral solution by algorithm GRA-LPi by , then GAP-LP1 is . Similarly, Gap-LP3 columns show the gap between the cost of our obtained integral solution and the optimal value of LP3. That is, GAP-LP3 is . Since all the LP relaxations are lower bounds of the optimum of the original problem, these two columns are upper bounds of the performance gaps between the obtained solution and the true optimal allocation. Note that when the problem size is large, e.g., and , we are not able to solve LP1, we use N/A to denote such cases.
In Table 1, we can see that there are several features of our proposed algorithm (GRA-LP3). First, although solving LP3 is not as efficient as solving LP2, it is still mostly tractable while the solution time of LP1 increases very fast and soon becomes intractable. On the other hand, the solution provided by GRA-LP3 could provide significant improvement over the solution that is obtained by using LP2. In the tests we presented, there are cases in which GRA-LP3 could produce the exact optimal solution, while there are only if one uses GRA-LP2. Therefore, we can conclude that GRA-LP3 could deliver higher-quality solutions than GRA-LP2 within reasonable amount of time.
Next we test our algorithm using a benchmark problem set AP (Australia Post), which was collected from a real postal delivery network in Australia, see . In  and , Ernst and Krishnamoorthy solve the p-hub location problems for the AP data set, and we test our algorithms using the hubs their solutions specified. In particular, some of the hub-to-city cost coefficients are non-symmetric in the AP data set. In our experiment, we make adjustment to it accordingly by specifying in-flow and out-flow coefficients separately for each . The results of our tests are shown in Table 2.
In Table 2, we test 15 AP benchmark problems. Since solving FHSAP-LP1 already produces optimal integral assignments for all 15 problems in less than 120 seconds, we omit it in the table. In our results, GRA-LP3 obtains optimal assignments for 14 out of the 15 test problems, and the cost is only higher than the optimal cost for the remaining one. Meanwhile. the time it spends to compute the solution is much less time than that of GRA-LP1. Therefore, we conclude that our approach is quite reliable and efficient in solving real problems.
5 Robust FHSAP
In previous sections, we studied the fixed-hub single allocation problem with deterministic demand. In practice, the decision maker may not have accurate demand information. In such cases, it is of great interest for the decision maker to have a robust policy, which protects him from any realization of demand. In this section, we propose a robust formulation for the FHSAP and provide an efficient algorithm for it.
We adopt the notations used in Section 2. However, instead of knowing the pairwise demand exactly, we only know that they are in the following set:
Here is the nominal demand, is a weight matrix and is the -norm () of a vector defined by
The right hand side in (5) is the “budget” of robustness, indicating one’s uncertainty level about the input data. Such an uncertainty set is quite common in the robust optimization literature with most common choices of to be , or . For a comprehensive review of the robust optimization literature, we refer the readers to Ben-Tal et al. .
Now we consider the robust formulation for FHSAP. We start from FHSAP-QP. In the robust formulation, we aim to minimize the worst-case cost for any demand realization that is in the set . Therefore, the robust formulation for FHSAP-QP can be written as:
One feature of this robust formulation is that given a set of , the inside maximization problem has an explicit optimal solution. Define
the inside problem can be written as
By using standard Lagrangian method, one can obtain the optimal value to (8) as:
where . Therefore, the robust counterpart of (8) can be written as:
Now if the inter-hub costs are the same (w.l.o.g., all equals to one), we can write as . Further relaxing the binary constraints on , we obtain a convex optimization problem:
In general, the optimal solution in (5) is fractional. However, as we will show in the simulation results, in most of our test problems, an integral solution is automatically obtained. And for the cases in which a fractional solution exists, we can again apply the geometric rounding technique introduced in Section 3 to obtain an integral solution.
5.1 Numerical Experiments
In this section, we perform numerical tests using the robust approach we proposed above. We show that the robust approach can indeed guard against data uncertainty without compromising much the average performance (both in terms of computational efficiency and solution quality).
In the following we use the robust approach (7) with . In such cases, the robust counterpart (5) will be a second-order cone program (SOCP). In the tests, we study different experimental setups. For each setup, we consider the following two approaches:
Simply use the nominal demand to obtain the solution. We use the GRA-LP3 approach discussed in Section 3. The integral solution obtained (after rounding) is denoted by .
We use the robust approach with a predetermined budget of robustness . We then solve the program (5) to obtain the optimal solution (combined with the same geometric rounding procedure if the solution is fractional). We denote this solution by .
To evaluate the solutions, we use two performance measures. First, we evaluate the solutions at the nominal demand. We write to denote the cost of using at the nominal demand and define as to measure the percentage gap between the costs in nominal cases. Second, we compute the worst-case costs of our solutions for the demand in the robust set we defined. That is, we compute the objective values of (5). We use to denote the worst-case cost of using and define as .
In the experiment, we test the 15 AP benchmark problems. We use the nominal demand of AP data set as the parameter which could be interpreted as the mean of the uncertain demand. The parameter is generated from a standard lognormal distribution, and then multiplied by . The parameter varies with the magnitude of demands in different benchmark cases. The results are shown in Table 3.
In Table 3, Time1 is the computational time of the regular approach and Time2 is the computational time of the it robust approach. First we can observe that although Time2 is larger that Time1, it is still tractable. Regarding the performances of the two approaches, we can see that Gap1 is smaller than in 11 cases out of 15 and is less than in 12 cases while Gap2 is larger than in 9 cases in 15 and is almost in the -city--hub case. Therefore, the additional cost for the solution that one has to pay at the nominal demand is much smaller than the potential additional costs one needs to pay if one use but the demands turn out to be adverse. Or in other words, the benefit of the robust approach outweighs the cost.
In this paper, we studied the fixed-hub single allocation problem. We made two contributions. First, we proposed a new solution approach for this problem by establishing a new LP relaxation formulation. This new LP relaxation lies in between two known relaxations in the literature, and by showing numerical results of this relaxation, we show that it has a good balance between computational complexity and the solution quality. Second, we propose a robust version of the FHSAP problem. The robust problem aims to minimize the worst-case cost when the demand is only known to be in a certain set. We propose an algorithm to solve the robust FHSAP problem and show that indeed the robust formulation can guard against data uncertainty with relatively little cost. We believe that both of our contributions may help people to find the desired model/approach when facing such problems in practice.
-  A. Ben-Tal, L. El Ghaoui, and A. Nemirovski. Robust Optimization. Princeton Press, 2009.
-  J. Campbell. Integer programming formulation of discrete hub location problems. European Journal of Operational Research, 72(2):387–405, 1994.
-  J. Campbell. Hub location and the p-hub median problem. Operations Research, 44(6):925–935, 1996.
-  J. Campbell, A. Ernst, and M. Krishnamoorthy. Hub location problems. In Z. Drezner and H.W. Hamcher, editors, Facility Location: Application and Theory. Springer, 2002.
-  C. Chekuri, S. Khanna, J. Naor, and L. Zosin. Approximation algorithms for the metric labeling problem via a new linear programming formulation. In SODA’01: The twelfth annual ACM-SIAM symposium on Discrete algorithms, pages 109–118, 2001.
-  A. Ernst and M. Krishnamoorthy. Efficient algorithms for the uncapacitated single allocation p-hub median problem. Location Science, 4(3):139–154, 1996.
-  A. Ernst and M. Krishnamoorthy. Solution algorithms for the capacitated single allocation hub location problem. Annals of Operations Research, 86:141–159, 1999.
-  D. Ge, S. He, Y. Ye, and J. Zhang. Geometric rounding: Dependent randomized rounding scheme. Journal of Combinatorial Optimization, 22(4):699–725, 2011.
-  J. Kleinberg and E. Tardos. Approximation algorithms for classification problems with pairwise relationships: Metric labeling and markov random fields. Journal of the ACM, 49(5):616–639, 2002.
-  J. Klincewicz. Heuristics for the p-hub location problem. European Journal of Operational Research, 53(1):25–37, 1991.
-  M. O’Kelly. A quadratic integer program for the location of interacting hub facilities. European Journal of Operational Research, 33(3):393–402, 1987.
-  M. O’Kelly, D. Bryan, D. Skorin-Kapov, and J. Skorin-Kapov. Hub network design with single and multiple allocation: a computational study. Location Science, 4(3):125–138, 1996.
-  M. O’Kelly, D. Skorin-Kapov, and J. Skorin-Kapov. Lower bounds for the hub location problem. Management Science, 41(4):713–721, 1995.
-  D. Skorin-Kapov, J. Skorin-Kapov, and M. O’kelly. Tight linear programming relaxations of uncapacitated -hub median problems. European Journal of Operational Research, 94(3):582–593, 1996.
-  J. Sohn and S. Park. The single-allocation problem in the interacting three-hub network. Networks, 35(1):17–25, 2000.