Abstract

We consider a capacitated hub location-routing problem (HLRP) which combines the hub location problem and multihub vehicle routing decisions. The HLRP not only determines the locations of the capacitated -hubs within a set of potential hubs but also deals with the routes of the vehicles to meet the demands of customers. This problem is formulated as a 0-1 mixed integer programming model with the objective of the minimum total cost including routing cost, fixed hub cost, and fixed vehicle cost. As the HLRP has impractically demanding for the large sized problems, we develop a solution method based on the endosymbiotic evolutionary algorithm (EEA) which solves hub location and vehicle routing problem simultaneously. The performance of the proposed algorithm is examined through a comparative study. The experimental results show that the proposed EEA can be a viable solution method for the supply chain network planning.

1. Introduction

In supply chain management, the design of distribution networks is one of the most important problems because it offers a great potential to reduce costs and to improve service quality. An important aspect of designing a distribution network is the determination of the locations of facilities such as warehouses, depots, and distribution centers. Classical facility location models assume that each customer is served on a straight-and-back basis on a given route while computing distribution cost. This situation is true only if the demand of each customer is a full truckload. However, in many applications arising in practice, the demand of each customer may be less than a truckload such that multiple customers are served in a single route and distribution cost depends on the sequence of customers on the route [1]. In this case, to reflect accurately the distribution cost of routes within a location model, the location-routing problem (LRP) should be solved simultaneously.

The LRP deals with determining the location of facilities and the routes of the vehicles for serving the customers under some constraints such as facility and vehicle capacities and route length to satisfy demands of all customers and to minimize the total cost including routing costs, vehicle fixed costs, facility fixed, and operating costs. Location-routing models are especially necessary for systems where the time horizon for the facility location decisions is not too long, and location costs are comparable to the routing costs. For an extensive review and classification of LRP, see Nagy and Salhi [2].

This study is motivated by the observation made in the postal service. Hub and spoke structure is a well-known configuration implemented in postal systems. It provides economies of scale by consolidating the traffic flows at the hubs as opposed to connect directly each origin-destination (O-D) point. In the standard LRP, we assume a structure of facilities serving a number of customers, who are connected to their depot by means of vehicle tours. No routes connect facilities to each other. On the other hand, this study deals more complex network structure than the standard LRP. We consider a variant of the LRP known as the hub location-routing problem (HLRP). In the HLRP, flow of packages, mail, or passengers from different nodes are collected at hubs, transferred between hubs along hub links in order to economically consolidate flows on the route, and distributed to their destinations. In a typical hub and spoke network, the flow of each O-D pair consists of three components, in other words, collection, transfer, and distribution, respectively. For handling the HLRP it is convenient to distinguish between two types of decisions, that is, hub location problem (HLP) and multihub vehicle routing problem (MH-VRP). The HLP is concerned with determining which nodes in a network are designated as hubs and which nonhub nodes are allocated to each hub. In -hub median problem, the total number of hubs is fixed as and the locations of number of hubs are determined such that the sum of the flow costs is minimized. The HLP can be categorized into single allocation and multiple allocations hub problem. In the single allocation, each nonhub node must be allocated only one hub whereas, in the multiple allocations, each nonhub node may be allocated to more than one hub. Further, the HLP can be categorized as either capacitated or uncapacitated. The capacitated HLP considers capacity constraint on the amount of flow through hub. On the other hand, there are no capacity restrictions on the uncapacitated HLP. While the HLP that has received the most attention from researchers is the uncapacitated single allocation -hub median problem, this paper deals with the capacitated single allocation -hub median problem for the part of HLP. We prescribe that exactly hubs need to be chosen within a set of potential hub locations. Each candidate location of hub has capacity restriction and fixed cost for establishing it as hub. In addition, we determine the allocation of each customer node to exactly one hub among the determined number of hubs. The flow of each O-D customer pair should go through one hub and at most through two hubs. Once the HLP is solved, the HLRP reduces to the MH-VRP. A set of homogeneous vehicles with the same capacity are available. Each customer must be served by exactly one vehicle for pickup and delivery service, respectively. The loads of each route cannot exceed the vehicle capacity and the total demands of customers allocated to each hub also should not exceed the capacity of that hub. We determine the routes of vehicles so as to minimize the routing costs and fixed vehicle costs.

Kuby and Gray [3] consider the hub network design problem where they determine the least cost set of direct or stopover routes for the traffic from a given set of points to a predetermined hub. Aykin [4] addresses the HLRP in continuous solution space where the hubs can be located anywhere in the service region encountered in air cargo transportation. He presents an algorithm solving a series of shortest path problems and a multifacility location problem in an iterative manner. Even though these two studies consider routing possibilities, their approaches are restricted in parallel with their assumptions. Nagy and Salhi [5] introduce the many-to-many location-routing problem similar to the HLRP where several customers send goods to others and the locations of each hub are to be determined. In their study, a hierarchical heuristic solution framework based on the concept of nested methods is presented. Bruns et al. [6] consider a problem arising in the parcel delivery operations of a postal service. They determine the locations of delivery hubs and the allocation of customer area to delivery hubs. In the case of routing costs, they present a route length estimation formula. The problem studied by Wasner and [7] is closely related to the HLRP. However, all interhub flow must go through a central hub rather than allowing all hubs to be directly connected to each other. Çetiner et al. [8] consider the combined hubbing and routing problem in postal delivery systems and develop an iterative two-state solution procedure for the problem in order to produce a route-compatible hub configuration. Catanzaro et al. [9] investigate the partitioning-HLRP, a peculiar version of the HLP involving graph partitioning and routing features and mainly arising from the deployment of an internet routing protocol. The partitioning-HLRP consists of partitioning a given network into subnetworks, locating at least one hub in each subnetwork, and routing the traffic within the network at minimum cost.

The HLRP is very difficult to be solved since it is composed of two NP-hard problems, the HLP and the MH-VRP. Therefore, the relatively large number of papers devoted to develop the approximation methods and the heuristic approaches can be classified into three categories such as clustering-based, iterative, and hierarchical heuristics, respectively. Clustering-based methods begin by partitioning the customer set into clusters, one cluster per potential depot or one per vehicle route. Then, they locate a depot in each cluster and then solve a VRP for each cluster [1013]. Iterative heuristics decompose the problem into its two subproblems. Then, the methods iteratively solve the subproblems, feeding information from one phase to the other [1417]. Hierarchical heuristics considers the location problem as the main problem and routing as a subordinate problem where main algorithm is devoted to solving the location problem and refers in each step to a subroutine that solves the routing problem [18, 19]. To solve the problem, the HLRP is formulated as a 0-1 mixed integer programming model. The integer programming model can provide the optimal solutions for small sized problem instances, but it is not practical for problems of large size in a reasonable computational time. Therefore, a novel solution framework based on an endosymbiotic evolutionary algorithm (EEA) is proposed which solves the HLRP and the MHVRP simultaneously.

The remainder of this paper is organized as follows. In Section 2, a mathematical model of the HLRP is developed. Section 3 presents an EEA and the detailed descriptions. Section 4 provides the computational results. Finally, we present the conclusion and the discussion of some future research directions.

2. Mathematical Model

The HLRP under study can be formulated as a 0-1 mixed integer programming model which is an extension of the many-to-many proposed by Nagy and Salhi [5]. The HLRP has some different characteristics with the many-to-many LRP of which exactly hubs have to be chosen within a set of potential hub locations and each candidate hub has a capacity restriction. In addition, the vehicle route of pickup and delivery purpose for each hub can be constructed separately as shown in Figure 1. Then the HLRP under study can be formulated mathematically using the following notations:: the set of nodes, ,: the set of candidate hub locations,: the set of customer nodes,: the set of type of vehicle route = ,: the fixed cost of establishing candidate hub node as a hub, ,: the capacity of loads associated with candidate hub node , ,: the distance between nodes and , ,: the amount of flow between origin and destination , ,: the number of hubs to be opened,: the transfer coefficient between two hub nodes,: the routing coefficient among the customer nodes,: the fixed vehicle operating cost,: the maximum capacity of each vehicle,: , if a route type can be constructed directly from node to node ; ,: , if candidate hub is used as a hub; ,: , if customer is assigned to hub ; ,: , if the flow from origin to destination routed via hub and ; ,: the amount of flow from hub to hub ,: pickup load on vehicle just after having serviced customer ,: delivery load on vehicle just before having serviced customer .

Consider

The objective function of the given model is to minimize the total cost including flow costs between all the O-D pairs, hub, and vehicle fixed costs. Constraint (2) ensures that any vehicle leaving from any customer must visit exactly one customer as the next node, and constraint (3) guarantees that the number of entering and leaving vehicles to each node is equal. Constraint (4) ensures that each customer must be assigned to exactly one hub. Constraints (5)–(7) forbid the illegal routes, in other words, the routes which do not start and end at the same hub. Constraint (8) requires that each customer can be assigned to the hub when it is determined to be opened. Constraint (9) requires that the total delivery and pickup loads on any hub do not exceed the corresponding hub capacity. Constraint (10) and (11) guarantee that all the traffic between each O-D pair must be routed via the hub subnetwork. Constraint (12) requires that the total amounts sending from any hub to all other hubs are equal to the total sum of the sending amounts from all customers assigned to the hub. Similarly constraint (13) also requires that the total amounts receiving to any hub from all other hubs are equal to the sum of the receiving amounts to all customers assigned to the hub. Constraint (14) ensures that the amounts of flow from hub to hub are equal to the sum of the amounts sending from all customers whose transportation route includes hub and . Constraint (15) requires that the total number of hubs is . Constraints (16)–(21) are valid inequalities proposed by Karaoglan et al. [1]. Valid inequalities are one of the most practical ways to strengthen the linear relaxation of the formulations. These inequalities eliminate some practical solutions from the solution space such that a stronger lower bound can be obtained for the problem. Constraints (16) and (17) prevent generating the subtour among the customers assigned to the same hub. Constraints (22)–(25) are known as integrality constraints. The last two of constraints (26) and (27) are nonnegativity constraints.

3. Endosymbiotic Evolutionary Algorithm for HLRP

The problem under study is hard to be solved exactly when the number of nodes increases. Therefore, we propose an EEA based solution method to get the near optimal solution of the HLRP within a reasonable time. Kim et al. [20] introduce the concepts and underlying fundamentals of the EEA. Nowadays, the EEA based algorithms have been applied and proven to be very effective on combinatorial optimization problems. Kim et al. [21] apply an EEA to the balancing and sequencing problems in mixed-model U-shaped lines. Rocha et al. [22] propose a transgenetic algorithm based on the process of endosymbiosis for the biobjective minimum spanning tree problem. Yadav et al. [23] present a robust optimization technique based on EEA for a multistage, multiperiod logistics system. Shin et al. [24] provide a multiobjective symbiotic evolutionary algorithm to solve the hub and spoke location problem under multiple scenarios. Sun [25] develops an EEA for the capacitated asymmetric allocation -hub median problem. He shows that the EEA is a very valuable approach to design the postal delivery network. In the EEA, the original problem is split into subproblems, and each subproblem has a population consisting of a set of elements. An element in the population for a subproblem can be a part of a complete solution to the original problem. The selection of hub locations and the vehicle routing decision are regarded as subproblems in this paper. The EEA well manages two subproblems and efficiently acquires high quality solutions of the entire problem. The EEA is a type of symbiotic evolutionary algorithm that can consider multiple subproblems at the same time [26, 27]. When the original problem is interwoven by multiple subproblems and we want to solve the original problem as a whole instead of considering them individually, the EEA can be a good option to search for the solutions of multiple subproblems concurrently. The EEA considers each subproblem along with its symbiotic partners, which are corresponding subproblems, in the original problem during the evolution. The concurrent consideration of subproblems may reduce the probability that the algorithm dwells on the local optimum.

The EEA for the HLRP under consideration breaks down the original problem into two different subproblems and it also considers the original problem itself as combined problem of two subproblems; one is the subproblem for hub location (SP-H), another is the subproblem for vehicle routing (SP-V), and the last is the combined problem considering both hub location and vehicle routing (CP-HV). Two subproblems and the combined problems maintain their own population. The populations for the subproblems regarding hub location and vehicle routing in the proposed algorithm are referred as POP-H and POP-V, respectively. These two populations play the roles of corresponding symbionts in the endosymbiotic theory, and individuals in those populations represent partial solutions of the original problem. Both populations evolve in such a direction that corresponding symbionts from both populations cooperate with each other to find the better solutions to the original problem. Since individuals in POP-H and POP-V are merely partial solutions, only when corresponding symbionts are combined appropriately into the endosymbionts, evaluation of the solutions to the original problem becomes possible. The population for the combined problem (POP-HV) plays the role of the endosymbionts in this algorithm. An endosymbiont carries the genes of all symbionts, which are partial solutions. Therefore, individuals in POP-HV represent solutions to the original problem. The individuals in POP-HV compete for survival with new offsprings that are generated through the mating of individuals from POP-H and POP-V. Eventually, POP-HV evolves toward a better population that contains better individuals in terms of the original problem. Figure 2 shows the concept of the EEA for the HLRP under consideration and depicts the overall evolution procedure of the proposed EEA with the interactions among POP-H, POP-V, and POP-HV. In Figure 2, evolution means that the individuals in POP-H and POP-V are separated but cooperate with each other while evolving and the individuals in POP-HV compete with those in POP-H and POP-V, and at the same time they evolve by themselves while keeping up their combined form within the population. The detailed procedures of evolution are described in the following Section 3.1.

3.1. Overall Procedures of EEA

In the algorithm, individuals in POP-H and POP-V are separately working as partial solutions but cooperate with each other. To make a complete solution of the original problem, a pair of individuals from each subproblem is required. On the other hand, individuals in POP-HV represent solutions of the original problem and they compete with those created by combination of individuals from POP-H and POP-V. In general, the search efficiency of evolutionary algorithms is often influenced by population diversity. An evolutionary search with a diverse population can continue to produce new generations and thus avoid becoming trapped at local optima [26]. Therefore, maintaining diverse populations is necessary for the long-term success of any evolutionary algorithm. Kim et al. [20] suggest a two-dimensional structure of toroidal grid as a population topology and the localized interactions among the populations instead of randomized and scattered interactions in order to maintain population diversity. They show the effectiveness of this policy in their previous work. For appropriate localized interactions in the processes among populations, we adopt the same population structure and the topological locations of corresponding individuals are defined as follows.

Each population forms a two-dimensional structure of toroidal grid with the same number of individuals, , and individuals in the population are mapped into the cells of the grid. Let , , be the location index in the toroidal grids. An individual in the population has its own location index and is surrounded by 8 neighbor individuals. The individual has corresponding individuals in the same geographical location at the toroidal grids of other two populations. Hence, when an arbitrary location is selected at a generation, the neighborhoods of individuals including at center of the 3 × 3 grid and 8 neighbor individuals in the POP-H, POP-V, and POP-HV are defined as , , and , respectively. Only the individuals in these three sets of neighborhoods are considered for the interactions among the three populations at the generation. The neighborhoods of from POP-H and from POP-V cooperate to find a good solution of the problem. Since each neighborhood contains 9 individuals, 81 (9 × 9) combinations are considerable as the candidate solutions of the original problem. The best combination among them is compared with the current best solution in the algorithm that competes with 9 individuals in from POP-HV as well. Based on the interactions among the sets of neighborhoods, the parallel search with partial solutions from POP-H and POP-V in subproblems and the integrated search with complete solutions from POP-HV in the combined problem are carried out simultaneously in all generations. The procedures of the proposed algorithm proceed as follows.

Step 1 (initialization of population). Generate initial individuals of POP-H, POP-V, and POP-HV, randomly. Set the best solution value, .

Step 2 (construction of neighbors). Select an arbitrary location and set up the neighborhoods, , , and , in each population.

Step 3 (cooperation between subproblems). Consider the following:Step : evaluate the fitness of all possible combinations that can be produced by the concatenation of individuals in and .Step : let be the best combination among those evaluated in Step , which becomes a candidate endosymbiont.Step : if , then update and keep as the best solution from subproblems for the competition in Step .

Step 4 (competition between entire problem and subproblems). Consider the following:Step : evaluate the fitness of individuals in .Step : let be the best individual and let be the worst one in .Step : if , then update and keep as the current best solution.Step : if , then replace with in .

Step 5 (evolution). Consider the following:Step : perform the evolution of individuals in by genetic_evolution (, , , ).Step : perform the evolution of individuals in by genetic_evolution (, , , ).Step : perform the evolution of individuals in by genetic_evolution (, , , and ).Step : release the evolved neighbors, , , and , to the population, POP-H, POP-V, and POP-HV, respectively.

Step 6 (iteration and termination). If the termination criteria of evolution are met, then stop. Otherwise, go back to Step and repeat the process.

In Step , not only POP-HV but also POP-H and POP-V evolve to find better solutions through competition and cooperation in the algorithm. Only individuals in the neighborhood in three populations are considered to perform the evolution. The procedure of genetic_evolution used in Step is as follows.

Genetic_Evolution (, , , and )

Step 1 (selection and crossover). First, select individuals from P using a selection mechanism and a random number is generated for each selected individual. When the random number is smaller than or equal to (crossover rate), the individual becomes a candidate of crossover operation. Then, create offspring by applying a crossover operator to a pair of candidate individuals.

Step 2 (creation of new generation). Individuals for replacement are selected from P and the selected individuals are replaced by the newly created offspring.

Step 3 (mutation). Select (mutation rate) individuals from P and apply mutation operator to each individual.

Step 4 (evaluation and update). Evaluate the newly generated offspring and update when the fitness of the best individual is higher than .

Each population requires its unique genetic representation for the corresponding problem and goes through its own genetic operations to evolve. The details of genetic representations and operators are described in the following sections.

3.2. Genetic Representations

The solution representation method is very critical to the successful implementation of EEA to the application area. The number of hubs which will be open must be specified as at the beginning of the search process. The proposed EEA contains two different subproblems, SP-H and SP-V. Each subproblem requires different genetic representation method due to different characteristics. In this paper, one-dimensional array is used to genetically represent individuals for SP-H. According to the characteristics of SP-H, the length of an individual is equal to the number of the candidate hub locations plus the number of customer nodes. Assume that there are three potential hub locations numbered from 1 to 3 and is given. Further, the locations 1 and 3 are selected as hubs and there are ten customer nodes in the network. Then, binary representation for hub selection and permutation representation for node allocation are applied to this study as shown in Figure 3.

After the hub locations are determined, each customer node is successively allocated to the hub. For the genetic representation of node allocation, problem specific encoding/decoding method where an individual is encoded as a permutation of node numbers and then decoded as one-dimensional array of hub numbers is proposed. In Figure 3, an individual for node allocation can be represented as (4 7 1 5 8 10 2 3 6 9) which is not an allocation yet. Then each customer node is considered one by one and assigned to the hub which has the maximum value of , defined as (28). Intuitively, each customer node is likely to be allocated to the hub which incurs the low vehicle routing costs with the already allocated customer nodes:where is the set of customer nodes which are allocated to hub in the th individual. Note that also includes the hub itself. This procedure continues until all customer nodes are allocated to the hubs. The proposed problem specific decoding method guarantees that each decoded solution does not violate the capacity constraints on hubs.

The MH-VRP considered in this paper assumes that all the vehicles have an identical capacity. All customers can be visited only once by a single vehicle for pickup and delivery respectively. When a new route of vehicle is created a fixed cost per vehicle is considered in order to minimize the total number of vehicles used. In order to solve the MH-VRP, each route of the vehicle should be complied with the clusters obtained in the corresponding individual of POP-H. Therefore, a permutation of customer node numbers is used to genetically represent individuals for MH-VRP. And then each individual is rearranged according to the clusters of corresponding symbiotic partner. Let be the number of customers who are assigned to hub . The following equation illustrates the genetic representation of an individual for SP-V:

The genetic representation of each cluster is encoded as a permutation of customer nodes already assigned to the hub. Its route is constructed by incrementally selecting customers in consideration of the demand of each customer and the loading capacity of vehicles until all customers have been visited. Initially, each vehicle starts at each hub and the set of customers included in its tour is empty. The unassigned customers are considered one by one in the order they appear in the permutation and the storage capacity of the vehicle is updated before another customer is selected. The vehicle returns to the depot when the capacity constraint of the vehicle is met or when all customers are visited. This decoding method is applied to each cluster successively. It ensures that the feasibility conditions are not violated by the crossover operations such that every permutation results in a feasible solution.

The proposed EEA consists of not only subproblems (SP-H and SP-V) but also a combined problem (CP-HV). An individual in POP-HV is a candidate solution for the original problem and consists of two different parts of genomes; one has an identical form of the individual in POP-H and the other has an identical form of the individual in POP-V. Therefore, an individual in POP-HV is a completely formed solution for the HLRP under consideration. The genomes for the hub location and vehicle routing are concatenated in a row. The genome for the hub location is placed at the head of the individual and is followed by the genome for the vehicle routing. The genome at the head and the genome at the tail use identical genetic representations for SP-H and SP-V which are already explained in this section. The concatenation of two different genes results in a one-dimensional array, but two separate genes within an individual in POP-HV are heterogeneous. Therefore, the genetic operations for POP-HV are carried out separately in two steps: one for the head and the other for the tail. Figure 4 shows the genetic representation of individuals in POP-HV.

3.3. Genetic Operators

At the end of each generation of EEA, only the individuals in the three neighborhoods, , , and , are considered for genetic operations such as selection, crossover, and mutation to generate new offspring for the next generation. Since a pair of individuals from each subproblem is needed to form a complete solution for the original problem, the individual in the pair shares the fitness value in the proposed EEA. The higher an individual has fitness value, the higher it has the chances to participate in mating pool. In this paper, a stochastic remainder sampling without replacement method [28] is used for selection operator. After selection process is completed, members of the newly generated chromosomes are mated at random. Then, the crossover operator is applied to produce two offspring. In this study, we use the position based crossover [28] for HLP as the following.

Step 1. Two individuals (P1 and P2) in the current neighborhood are selected at random.

Step 2. Select a set of positions from one parent (P1) at random. Produce a protochild by copying the genes on these positions into the corresponding positions of it. Delete the genes which are already selected from the other parent (P2).

Step 3. Place the genes into the unfixed positions of the protochild from left to right according to the order of the sequence in P2 to produce one offspring. This process is applied to each genome for hub selection and node allocation successively.

Further, we present a modified order crossover (MOX) operator for MH-VRP and its procedure can be explained as the following.

Step 1. Two individuals (P1 and P2) in the current neighborhood are selected at random.

Step 2. A vehicle number used in one parent (P1) is randomly chosen. The corresponding genes, which are the sequence of customer indices for the chosen vehicle, are copied to the beginning of the new offspring. The indices of already inherited genes from P1 are deleted from the other parent (P2).

Step 3. The remaining genes in P2 are copied into the offspring in the order in which they appear in P2. This process is applied to each cluster successively.

The proposed MOX operator is illustrated with an example as follows. Assume that there are 12 customer nodes of the first cluster only for simplicity, nodes 1, 5, 7, and 9 are assigned to vehicle 1, nodes 2, 4, 8, 11, and 10 are assigned to vehicle 2, and nodes 3, 6, and 12 are assigned to vehicle 3, respectively. Assume that the randomly selected vehicle in the parent 1 is vehicle number 1. Then possible generation of offspring can be described as shown in Figure 5.

In addition, a swap mutation operator is used which selects randomly two nodes in different clusters and swaps their assigned hubs. If the feasibility is violated by the mutation operator a penalty is imposed to the fitness function when the individual is evaluated.

To improve the quality of solutions for the various combinatorial optimization problems, various local search methods have been applied to the evolutionary algorithm. In this paper we adopt the 2-opt exchange method to improve the quality of solution for the MH-VRP. For each pair of customers in the route of each vehicle, the 2-opt exchange heuristic swaps their locations in the route and finds the best combination resulting the most great improvements. It has been applied to only the elite individual of the neighborhoods, , , and in each generation after crossover and mutation operators are finished.

4. Computational Results

This section evaluates the effectiveness of the proposed EEA through a comparative study. For this purpose, the proposed algorithm is programmed in C++ and run on a personal computer with a 3.6 GHz Pentium 4 CPU and 8.0 GB RAM. The parameters in EEA have been determined via the preliminary experiments. A 10 × 10 toroidal grid is used for POP-H, POP-V, and POP-HV. Thus the size of each population is 100. The crossover rate and the mutation rate are set to 0.8 and 0.05, respectively. To determine the crossover rate (p_crossover) and mutation rate (p_mutation), an experimental design is developed at different levels. We tested p_crossover values in the range of [0.5, 0.9] with an increment of 0.1 and p_mutation values of 0.01, 0.05, 1.0, 1.05, and 2.0. Finally we selected the pair, p_crossover = 0.8 and p_mutation = 0.05, which gives the highest performance. When the current number of generations reaches at a specified maximum number depending on the size of solution space, the proposed algorithm terminates.

Computational experiments are carried out to evaluate the effectiveness of the proposed EEA over multiple ant colony optimization (ACO) algorithm which is one of the most recent LRP heuristic by Ting and Chen [29]. We generated a test set that consists of a wide variety of problems. In the experiment, the number of customers was set at levels 100 and 200, whereas the number of candidate hubs was set at 10 and 20. Spatial distribution of customers may also affect the performance of the algorithms. It is considered by two factors: number of clusters and clustering ratio. The number of clusters is set at levels 0, 3, and 5, where level 0 refers to uniformly distributed customers. The clustering ratio is defined as the ratio of number of customers that belong to a cluster to the total number of customers. This factor is set at levels 50% and 100%. For all the instances, we generated customer demand from an uniform distribution in the range of . And vehicle capacity is set to 150 which corresponds to an average 10 customers on a route. The values of discount factors for transfer and routing are set to and . In all of these instances, hub fixed cost was set equal to 10.0 for all hubs and vehicle fixed cost to 1.0 for all vehicles. The comparison of computational results is given in Table 1.

In Table 1, columns 1, 2, 3, 4, and 5 indicate the problem number of test problems, the number of customers, the number of candidate hub locations, the number of clusters, and clustering ratio, respectively. The sixth and seventh columns show the solution values and CPU times (in seconds) obtained from the ACO algorithm. The results from the proposed EEA are reported in the eighth to tenth columns. The improvement ratio in column ten indicates the relative improvement ratio of solutions obtained by EEA over ACO. Each entry of the solution value obtained by EEA and ACO is the average of the 10 runs. Though we tried to solve the proposed 0-1 mixed integer programming model using Xpress-MP 7.6, commercial optimization software, we found that the mathematical formulation for the HLRP could not get the optimal solutions for the test problems with more than 10 nodes in a reasonable computational time. The mathematical programming approach is not practical to obtain a solution for the real sized problems.

The results of computational experiments indicate that the proposed EEA outperforms the existing ACO based solution method and generates improved solutions. As shown in Table 1, the proposed EEA provided better solutions in all test problems than the ACO algorithm. These encouraging results are due to the characteristics of EEA which it is originally designed to handle complex structures consisting of multiple subproblems. The distinguishing feature of the proposed EEA is that it maintains endosymbionts that are a combination of an individual and its symbiotic partner. The existence of endosymbionts can accelerate the speed that individuals converge to good solutions. It seems that the EEA can take advantage of a proper balance between the parallel search with partial solutions and the integrated search with entire solutions. This enhanced capability of exploitation together with the parallel search capability of traditional symbiotic algorithms results in finding better quality solutions than existing hierarchical approach based on ACO algorithm. As distribution costs account for a large portion of the company costs in most industries, it is readily justified to invest more computation time in order to improve the solution quality even by a small amount. It can be observed that the computation times required by both algorithms increase by approximately the same factor as the number of customers increases.

In addition, the progresses of the best solutions in both algorithms have been compared. For the experiment, problem number 16 in Table 1 is used. It consists of 200 customers and 20 candidate hub locations. Figure 6 shows the progress of the best solutions in each generation for the test problem from the EEA and ACO algorithm. Lines in the figure represent the change in the value of the best solutions in both approaches during the evolution.

5. Conclusion

In this paper we considered the capacitated HLRP which is composed of two subproblems, HLP and MH-VRP, encountered in postal delivery service. We developed a 0-1 mixed integer programming model and proposed a new integrated approach based on EEA in order to solve the combined problem simultaneously. In the HLRP, the hub selections and customer allocation decisions will affect the total cost of vehicle routes and the structure of vehicle routes will also affect the decisions of hub selections and allocation of customers to hubs. The proposed EEA iterates between location and routing phases in order to search for better solutions. The experimental results show that the proposed EEA outperforms the existing heuristic. The encouraging results suggest that the proposed EEA can be a valuable approach to design the postal delivery network and it can be applied to other combinatorial problems that contain multiple levels of decision making.

In the future research, a parallel implementation of the EEA can be considered for the rapidity of convergence in search process. Another future work is the development of robust optimization scheme in the environment of demand uncertainty.

Conflict of Interests

The author declares that there is no conflict of interests regarding the publication of this paper.

Acknowledgment

This work was supported by Hankuk University of Foreign Studies Research Fund.