Abstract
Hubs are critical facilities in the power projection network. Due to the uncertainty factors such as terrorism threats, severe weather, and natural disasters, hub facilities may be disrupted randomly, which could lead to excessive cost or loss in practice. One of the most effective ways to withstand and reduce the impact of disruptions is designing more resilient networks. In this paper, a stochastic programming model is employed for the hub location problem in the presence of random hub failures. A heuristic algorithm based on Monte Carlo method and tabu search is put forward to solve the model. The proposed approach is more general if there are numbers of hubs that would fail even with different failure probability. Compared with the benchmark model, the model which takes the factor of stochastic failure of hubs into account can give a more resilient power projection network.
1. Introduction
Power projection is a term used in military to refer to the capacity of a state to apply national transportation network to rapidly and effectively deploy and sustain forces in and from multiple dispersed locations to respond to crises, to contribute to deterrence, and to enhance regional stability [1]. In the process of power projection, persons and materials are consolidated at the transportation hubs (station, port, airport, etc.) through motorized march firstly, and then they are delivered using backbone transportation (highway, railway, marine transport, etc.). Finally they are deployed to dispersed locations through motorization march according to the mission demands [2]. The ability of a state to project its forces into an area may serve as an effective diplomatic lever, influencing the decision-making process and acting as a potential deterrent on other states’ behavior. Examples of power projection include the US-led Invasion of Iraq [3] and the Russians invasion of Georgia [4]. Another typical example of the power projection network is the earthquake relief action of the Chinese army in 2008 Sichuan earthquake. About 160 thousand troops and eight million tons of materials were project from thousands of kilometers to the disaster area by road, rail, and air in three days after the earthquake [5].
Generally, a hub-and-spoke structure is adopted in the power projection network. Hub facilities play key roles in the network by concentrating, distributing, and switching traffic flows instead of transferring flows between each origin–destination (O–D) pair directly [6]. There are two basic assignment rules for the connection of nonhub nodes and hub nodes, which are entitled by single-allocation (SA) and multiple-allocation (MA) rules [7]. The SA rule refers to each nonhub node to connect with only one hub so that all flows from an origin must travel to the same hub. In contrast, the MA rule is more flexible in routing by allowing each node to interact with more than one hub [8]. A typical real-life example of multiple-allocation is the express delivery system. Packages from one origin are assigned to different hubs depending on their destinations. The location of hubs has significant influence on the efficiency and safety of the network. Suffered from the uncertain risk of enemy’s attack, terrorism threats, severe weather, natural disasters, and so forth, hubs would be disrupted randomly and make the network partially or completely unavailable. Resilience is the ability of a network to withstand and reduce the impact of disruption events [9, 10]. Design of a resilient power projection network is a very important practical issue for strengthening the ability of strategic power projection. For response to hub failure, two strategies are usually adopted which include reactive (e.g., repairing [11] and rescheduling [12, 13]) and proactive strategies (e.g., adjusting the network hub location strategy [14, 15] and protecting key hubs [16–19]). It is sensible to take into account the impact of hub random disruptions in advance for designing a more resilient network.
A number of studies addressed hub location problem in the presence of hub disruptions. Bi et al. [2] assumed that the enemy would select the hubs for attack which result in the greatest damage to the network. They gave a resilient hub location strategy considering the failure of specific hubs. However, the enemy usually could not have perfect information about the location of hubs. Other risks such as severe weather, natural disasters, and traffic congestion are uncertain. It is more common that hubs fail randomly. An et al. [20] did some analytical work on the reliable hub-and-spoke design with consideration of the stochastic failures of hubs. They employed nonlinear mixed integer programming models and used Lagrangian relaxation and branch-bound method to solve the models. Azizi et al. [21] incorporated hub unavailability into the classical single-allocation p-hub median problem. They assumed that once a hub stopped normal operations, the entire demand initially served by this hub was handled by a backup facility. They proposed a mixed integer quadratic programming model and a metaheuristic algorithm. The latter two studies both assumed that there was at most one-hub failure in the network, and it was somehow unreasonable in practice. For example, if there is a 10-hub network with hub reliability 0.7, the expected number of hub failures is 3 rather than 1. In addition, the proposed approaches in the latter two papers which translated the failure probability directly into the mathematical expectation are not suitable for the situation where hubs have different failure probability. More relevant works can be found in the review papers of Reggiani et al. [22] and Mattsson and Jenelius [23].
Note that the reliable server assignment problem is similar to the resilient problem defined in this paper. However, there are some distinctions between resilience and reliability. Resilience refers to the capacity of networks to absorb and return to normal conditions after a critical event, while reliability refers to the probability of a device performing its purpose adequately for the period of time intended under the operating conditions encountered. For example, connectivity-based reliability measures assume that a network is not functional even if a single node fails or becomes disconnected. In practice, however, a network continues to serve the remaining connected components even though one or more nodes have become disconnected. Although network reliability measures have been frequently used in the literature as network design criteria, network resilience measures represent a viable alternative approach. The network resilience measure used in this work is more comprehensive than traditional network reliability measures since it incorporates traffic demand into its determination. For more discussions about resilience and reliability, see Mattsson and Jenelius [23].
In fact, when the hubs in the network fail randomly, it is hard to determine which one of the hubs would be disrupted in the future. Sometimes it is also out of control to find a suitable backup hub and alternate route. Therefore, compared with traditional hub location problems, it is more complex for resilient hub location in power projection network considering random hub failures. In this paper, a stochastic programming model is employed for the hub location problem in the presence of random hub failures. A heuristic algorithm based on Monte Carlo method and tabu search is put forward to solve the model. The superiority of the proposed approach is illustrated by a computational example compared with the benchmark model.
2. Problem Description and Assumptions
The main work in this paper is determining the location of hubs and the allocation relationship of nodes to hubs. First, we assume that the network is fully connected and all the nodes are in normal state. Second, we select some nodes as hubs randomly and determine the best assignment for nonhub nodes and hub nodes. Once the allocation relationship of nodes to hubs is determined, the redundant links are removed. Then we get a hub-spoke distribution network, as shown in Figure 1. Persons and materials are consolidated and deployed through the hubs. It is assumed that when hub failure occurs, its function will be lost and cannot be recovered limited time. However, the traffic generated at this hub and the flow through this hub will not be lost; in other words, this hub point changes into a nonhub point. The traffic flow will be replanned based on the shortest path principle through the remaining hubs. For example, assume that the traffic flow between is originally transported via the hubs . When hub fails and the other hubs are normal, the traffic flow can select route for two-hub stop transport or route for one-hub stop transport. The specific route to be used depends on the transport cost. Due to the fact that nonhub nodes cannot be directly connected, if a pair of start-end point cannot find a suitable hub transport, then the traffic loss occurs. The cost of loss can be measured by the cost of direct transport multiplied by a penalty factor.
Other basic assumptions are as follows:(1)Hubs are connected by a complete network.(2)A nonhub node can be connected with multiple hub nodes, but nonhub nodes cannot be directly connected.(3)A discount factor exists for the transportation cost when backbone transportation is used between the hubs.(4)The hubs are randomly disrupted, but the reliability of the hub (corresponding to the probability of failure) is known.(5)The hubs and links are uncapacitated.(6)The route planning and emergency scheduling follow the shortest path principle.
The resilience of the power projection network is defined as follows [2]:where is the resilience of the power projection network, is the transportation costs of the network in normal situation, and is the transportation costs of the network when some hub nodes have been disrupted.
3. Model Parameters and Variables
The model parameters and variables mentioned in this paper are defined as follows” : set of nodes : set of links : directed graph of the network : set of hub nodes : number of hubs : flow matrix, being the flow from to , and : unit transport cost from node to , and : cost discount factor between hubs : unit transportation cost from node i to j through hub pair : reliability of the hub , : reliability of the hub pair (if , , otherwise ) , 0-1 state vector for hub set (, hub normal; , hub failed) : a random variable, indicating the probability of finding the backup hub and alternate route when a node fails to pass through the hub, : penalty factor, penalty cost caused by goods transportation between nodes failed due to hub disabled ( in this paper) : decision variable, whether node is selected as a hub or not : decision variable, whether the transportation from node to go through the hub pair : the decision variable, whether the transportation from node to goes through the backup hub pair : Monte Carlo repetitions : total transportation cost for the th hub in a specific state given the probability of hub failure : the expected transportation cost estimation of projection network under a failure probability using Monte Carlo method
4. Models and Algorithm
In order to illustrate the importance of considering the possible failure of the hub in the design of the network, the multiple assignment hub location model in the normal state is provided as the benchmark model, and then a stochastic programming model considering the probability of hub failure is proposed.
4.1. Benchmark Model
The uncapacitated multiple assignments hub location model [24] for the power projection network without hub failure is as follows. The model is used to determine which nodes will be selected as hubs and the relationship of the nonhub nodes and hubs.
The objective function (2) minimizes the total cost of the power projection network in the normal state without hub failure. Constraint (3) requires that there are hubs to be located. Constraints (4) refer to the route choice of each O-D pair through one hub or a pair of hubs. Constraints (5) and (6) indicate that each O-D pair can only be assigned to an existing hub pair. Constraints (7) are binary requirements.
4.2. The Stochastic Programming Model under Random Hub Failures
Based on the hypothesis and notations mentioned above, the stochastic programming model for hub location in the power projection network considering the probability of hub failure is as follows:
where is a random variable, which refers to the probability that a backup hub and alternate route can be found in the network when some hubs are disrupted. The objective function (8) represents the expected total transportation cost. The first item on the right side of the equation is the expected transportation cost when all hubs are workable, the second is the expected transportation cost when the traffic flow is affected by the hub failure, and the third one is the penalty cost of traffic loss when a viable alternate route cannot be found in the network. Formula (10) is routing constraint, which guarantees that the transportation for any pair of nodes can only be executed on one route under a specific state. Constraints (11) and (12) indicate that the backup route selection only occurs when the pair of hubs is still workable during routes replanning. Formula (13) is 0-1 integer constraint.
The two above-mentioned models are interrelated. When the reliability of hub pair equals 1, that is, there is no hub failure, the last two items in (8) will be equal to 0. In this case, model (II) can be simplified to model (I).
4.3. Estimations of the Expected Projection Cost and Resilience Using MC Method
When the hub in the power projection network randomly failed with a certain probability, the Monte Carlo method [24] can be utilized to estimate the expected transportation cost based on the Floyd shortest path algorithm, and then the network resilience can also be calculated. The method is as follows.
Denote as the initial hub set in the power projection network, as the structure diagram, and , is the edge of . The length of is defined as follows:
Let the initial weight matrix of the network be and the routing matrix be and represents the subscript of the first hub in the shortest path from node i to node j. When there is no hub failure, the shortest path and unit transport cost matrix for any start-end point can be calculated based on the Floyd shortest path algorithm [25]. The specific steps are as follows.
Input weight matrix: , .
Calculate , and ; calculate , and , if .
The element in is the weight of the shortest path between any two points, the element in is the subscript of the head for the first arc of the shortest path from node i to node j, and the relationship between the nonhub and the hub can be obtained from the matrix. The shortest path for is .
The overall cost matrix under the optimal allocation route can be obtained by the above-mentioned iterative process when the hub is workable, and the total transportation cost of the network is calculated as
Considering the situation where the hub is randomly defeated with a certain probability, take as the 0-1 state vector for the hub set. When the hub is workable, ; otherwise, . A specific Monte Carlo experiment is designed and then repeated a total of times. During each repeat, p random numbers are generated using the uniform interval distribution. Denote these random numbers as and do comparison with the hub reliability . If , ; otherwise, . Let the state of the hub set be for the th experiment, and the normal working hub set be . After substituting into (14), the weight matrix of the power projection network at the th experiment will be obtained. Repeat the Floyd algorithm mentioned above, and calculate the transportation cost for the power projection network at the th experiment. The expected transportation cost of the power projection network with a failure probability obtained by the Monte Carlo method is calculated asAccording to the definition of resilience, the expected resilience of the power projection network with the specific probability of failure is calculated as
4.4. Tabu Search Algorithm
Model (II) is one of the combinatorial optimization problems with random variables. It is difficult to solve such problems by using traditional quantitative methods [7] and the programming software. In this paper we employ the tabu search [26] based on the Monte Carlo method to solve the model. The basic idea of the solution is to randomly generate a set of nodes, and then the former nodes are treated as the initial hub solution set. The neighborhood area for the solution can be constructed by using 2-swap exchange; that is, the swap of a nonhub and the hub can be used to construct the neighborhood area. Given the hub set, formula (16) for the expected network cost can be treated as the evaluation function.
In the process of optimization, the acceptance for one exchange of a hub node and a nonhub node in current solution is called a move. The moving rule is as follows: all solutions in the neighborhood area of the current solution are sorted according to the evaluation function and used as candidate solutions. If the move produced by the candidate solution with the smallest value of the evaluation function is not in the tabu table or in the tabu table but satisfies the aspiration criterion, it will be accepted and the current solution is updated for the next iteration; otherwise the next candidate solution gets the chance; if all candidate solutions are tabu, the best candidate solution is selected as the new current solution for the next iteration. It should be noted that the Monte Carlo experiment is carried out for each move.
The overall process of tabu search algorithm based on MC method is as follows.
Step 1. Set the tabu table to be empty, the maximum run times for algorithm to be Max Run, the maximum number of iterations for each run to be Max Iteration, the number of iterations for the optimal and consistent solution to be Max Count, and the node tabu frequency matrix to be the zero matrix.
Step 2. Generate the initial solution and take the initial solution as the current one and the current solution as the optimal one.
Step 3. Generate the neighborhood area for current solution. For each hub set in the neighborhood area, the corresponding evaluation function values are estimated by using the mentioned Monte Carlo experiment method and then sorted as the candidate solutions.
Step 4. Is the best candidate solution better than the optimal solution? If yes, turn to Step 7; otherwise, go to Step 5.
Step 5. Are all solutions tabu? If yes, select the best candidate solution as the current solution; otherwise, select the nontabu and the best candidate solution as the current solution and then turn to Step 6.
Step 6. Update the current solution, the tabu and tabu frequency, Iteration + 1, and Count + 1.
Step 7. Take the solution as the current one, and update the optimal solution, the tabu and tabu frequency, and Iteration + 1.
Step 8. If the number of iterations and the number of iterations for the optimal and consistent solution do not reach to the limitation, the iteration will continue. If the limit has been reached, the count of runs is incremented and go to Step 9.
Step 9. Eliminate the nodes with high frequency of tabu, and this results in a new initial feasible solution, and go to Step 2.
Step 10. Max Run is reached, and then the algorithm ends, and the result outputs.
5. Illustrated Example
Suppose that there are 20 projection nodes; the horizontal and vertical coordinates are sampled randomly from , as shown in Table 1. The spatial distribution of nodes is shown in Figure 2. The distance matrix between nodes is calculated from the Euclidean distance of the nodes on the plane, in kilometers. The traffic flow between nodes is generated by randomly getting value from in tones. The unit transportation cost is set at 0.6/RMB/ton/km. There are a total of hubs to be located, and the transportation cost discount coefficient of hub links is set to .
Model (I) and model (II) are solved by tabu search algorithm, and, based on Matlab 7.10, both run on Intel Core i5 2.20 GHz/8.00 GB computer. In tabu search algorithm, the tabu length is set to 15, the tabu frequency is set to 5, the maximum number of iterations is set to 50 for each run, and the number of iterations for the optimal and consistent solution is set to 20. The algorithm runs a total of 20 times. The parameter values in the tabu procedure are set using the experiences from our early related work; see Bi et al. [2]. The number of Monte Carlo trials for each hub reliability level is 10,000, and the statistical standard error is about 0.2%.
Model (I) is the conventional hub location scheme, and the model solves the result in 2.8 seconds. The optimal hub location scheme is 1 3 4 5 7 10 11 12 14 15. As shown in Figure 2 by symbol , the corresponding objective function value is RMB. Since model (I) is optimized without considering any hub failure, the optimal hub location strategy is the same corresponding to different hub reliability levels. When the power projection network is faced with uncertainty threat factors, the expected cost of the power projection network corresponding to Model (I) can be estimated by the Monte Carlo experiment method proposed in Section 4.3.
With the hub failures are taken into consideration in the early stage of the network design, the hub location scheme can be obtained by solving Model (II). Table 2 shows the hub location scheme, expected transportation cost, network resilience, and solution time cost for Model (II) under different hub reliability. The related results of Model (II) are compared with those of Model (I).
As shown in Table 2, when the hub reliability equals 1, that is, no hub fails, the hub location strategy in Model (II) is exactly the same as that in Model (I), indicating the correctness of Model (II). When the hub reliability is at a high level (), the hub location strategy corresponding to Model (II) has one-hub difference with the conventional hub location scheme. The expected cost and network resilience are almost the same, and Model (II) is a little better. This is due to the fact that the total expected transportation cost of Model (II) includes the expected transportation cost of the default route, the expected transportation cost of the emergency planning route, and the expected penalty cost of the traffic loss. When the hub reliability is high, the expected transportation cost of the default route will be the dominant part in the optimization process, and therefore the optimal hub location schemes provided by Model (I) and Model (II) are almost the same. With the decrease of the hub reliability, the probability of the hub random failure increases, and the difference between the hub set generated by Model (II) and Model (I) also increases. Meanwhile, the expected cost and difference in network resilience increase. The network designed with Model (II) is significantly superior to Model (I). This is due to the fact that when the number of failed hubs is more than expected, the proportion of the expected transportation cost for emergency planning route and traffic loss in the total expected cost increases, and then more effects will be made on the optimization process to get the optimal evaluation solution.
Figure 3 shows the comparison of network resilience between Model (II) and Model (I) at different hub reliability levels. It can be concluded that the network resilience of the hub location scheme designed by Model (II) is superior to the conventional hub location scheme provided by Model (I). With the probability of hub failure increase, the advantage is more obvious. Therefore, it can improve the resilience of the power projection network to a certain extent, when taking the factors of hub random failure and the hub reliability level into consideration in advance.
It should be emphasized that the solution time cost by Model (II) differs from ten minutes to more than three hours. The lower the hub reliability is, the longer the solution time costs. The reason is that solving Model (II) requires a large number of Monte Carlo experiments. Considering such great amount of calculations, the program can be improved in the optimization process. On one hand, since the nonhub can only be connected with the hub, the Floyd algorithm for a given hub set can just be iterated by p times. On the other hand, it can be seen from Table 2 that the difference in hub set is gradually increased with the decrease of hub reliability. So when using the tabu search algorithm, the optimal solution at a reliability level can be used as the initial solution for tabu search at next reliability level. These mentioned improvements can decrease the calculation burden greatly in the optimization process. In addition, the current computing speed and overall performance of computers are high enough, so that a large number of Monte Carlo experiments for parameter estimation become possible. The longest time for Model (II) to get the solution is 3.7 hours, and it is acceptable in the network design period.
6. Conclusion
In this paper, we study the hub location problem when random hub failure occurs in the power projection network. A stochastic programming model is employed account for the scenario. A heuristic algorithm based on Monte Carlo method and tabu search is put forward to solve the model. The superiority of the proposed model is illustrated by the change of the network resilience in response to the different probability of hub failure with the benchmark model. Results show that the factor of stochastic hub disruptions should be taken into account in advance for design of the power projection network, which can improve the resilience of the power projection network to deal with the hub random failure to a certain extent.
Compared with the situation that only some specific hub fails, the problem of hub random failure is more general. Dealing with the random hub failure, methods in the existing literature are only applicable to cases where there is at most one-hub failure. However, in this paper, the Monte Carlo method can be used to simulate the state of a hub set at any probability, which is applicable to the case where there are several hubs failure with different probability. Therefore, the approach proposed in this paper is more general.
Conflicts of Interest
The authors declare that there are no conflicts of interest regarding the publication of this article.