A Hybrid Algorithm for Location-Routing Sustainable Optimization under Fuzzy Demand

With the globalization of the supply chain and the change of demand environment, designing an effective logistic system in the sustainable development of the supply chain becomes more critical. This study proposes a location-routing problem to determine an efficient integration of single factory and multi-distribution centers and multi-customers in uncertain demands. This problem can be regarded as an optimization integrating location, distribution decision, and routing management. The objective function is to minimize the total cost and satisfy all the requirements, which is a highly complex NP-hard problem, so a hybrid algorithm of genetic algorithm (GA) and tabu search (TS) algorithm is proposed. A fuzzy c-means clustering algorithm is used to produce an initial solution. Fuzzy triangular number and confidence interval transformation are used to deal with fuzzy customer demand. The research findings concludes that (i) determine the numbers of facilities with locations that should be opened and (ii) minimize the total cost in supply chain. The experiments prove that the proposed hybrid algorithm of GA and TS algorithm overcomes the defect of local optimum in the literature viewpoint, and the optimization algorithms can effectively solve the location-routing problem.


I. INTRODUCTION
In order to remain competitive in the global manufacturing environment, effective supply chain strategic planning has become a vital factor. In integrated supply chain decision making, location problems (LP) and vehicle routing problems (VRP) are the two most studied areas in theory and practice. In LP, VRP, and other things, based on the flow decision model, location-routing Problem (LRP) arises. Therefore, transportation costs should be considered when making a facility location decision. By establishing the LRP model, for many in the case of customers and multiple facilities, it can simultaneously solve a series of related problem from determining the optimal number and capacity of distributions to seek the optimal transportation plan, route arrangement between overall issues, thereby reducing supply chain costs.

A. RELATED WORK
The location-routing problem determines the location of the depot from multiple candidate vehicle depots and finds the optimal vehicle schedule and route set based on the minimum traveling distance criteria [1]- [3]. In the 1980s, some scholars have studied heuristic algorithms for this problem. Van der Bruggen et al. [4] developed a simulated annealing algorithm, and Madsen et al. [5] proposed a parallel insertion's algorithm to settle time windows variable. Dumas et al. [6] used an accurate procedure to solve the multi-vehicle problem with time windows. The problem of power plant cycle positioning was proposed by Wan Jing et al. [7]. Osvald and Stirn (2008) [8] studied the distribution business of fresh food and applied the solution algorithm to the corresponding VRPTW, paying particular attention to the decline of product freshness. Similarly, a nonlinear mathematical model of production scheduling and vehicle routing considering time window is proposed by Huey Kuo Chen et al. (2009) [9]. Akkerman et al. (2010) [10] conducted a thorough survey of food distribution research, with a separate chapter devoted to transportation planning. Poorya Farahani et al. (2012) [11] developed an iterative solution, which integrated short-term production and location plans.
On this basis, some convergence schemes are proposed as interfaces for production scheduling and allocation problems. Nasab et al. (2013) [12] proposed a multi objective mathematical programming formula to explain this problem, the study used heuristic algorithm to a standard reference set VOLUME XX, 2017 1 to solve objective function and verified the effectiveness of the method. M. Shojafar (2015) [13] focused on better utilization of VM resources according to the capacity and performed fuzzy theory and GA to achieve better results. Y. Marinakis et al. (2016) [14] proposed iterative local search and variable neighborhood search algorithm under random demands to study the location routing problem. Because there are no benchmark instances for this form of problem in the literature, several have been the study created some new test instances based on capable location routing problem. Z. A. Akpunar et al. (2020) [15]  Previous studies focus on assuming that consumer demand obeys a specific probability distribution. However, past data on consumer demand may not always be reliable nowadays. Therefore, fuzzy demand better suits today's supply chain environment, and fuzzy set theory provides an alternative approach to deal with this uncertainty.

B. PAPER CONTRIBUTIONS
The proposed production/distribution problem under the fuzzy demand model attempts to minimize total costs. The dynamic location model studied in this paper takes into account a multi-period operating environment, with changing requirements from period to period. The problem can be solved in two stages. The first phase determines which distribution center (DC) will be enabled and which the enabled distribution centers will assign customers. The vehicle routing problem under fuzzy demand is solved in the second stage for an open distribution center. This problem is relatively simple to formulate with many solutions. Unfortunately, it is NP-hard. In this paper, fuzzy requirements are introduced, and hybrid algorithms are adopted to solve the problem.

C. OUTLINE
This study is organized as following: section 2 proposes the problem and formulates the original model. Subsequently, section 3 generates the initial solution. Next, section 4 presents a hybrid algorithm to solve the objective function. Section 5 gives an industrial case for implementing the feasibility of applying the proposed model. Finally, discussions and conclusions are drawn in section 6 and section 7.

II. MODEL OF LOCATION-ROUTING PROBLEM
In this paper, the production and location-routing problem has a plant, multiple distribution centers, and multiple customer positions, and the assumptions are following: • There is no stock in the factory; products are distributed directly to the distribution center after being produced. The production capacity is infinite in the factory.
• There are fixed start-up costs in the distribution center, which has capacity limits. Each distribution center can be supplied to multi-customers.
• Number of vehicles in each distribution center is known. Each car can only belong to one distribution center.
• Customer requirements vary within a certain range. In this fuzzy demand environment, the paper need to give an initial solution produced by density-based clustering and route solution. The cluster principle is the level of a customer node close to the distribution centers, and the optimized solutions can get more economical distribution route cost (cost of total distance). The first level operation mode is delivering the goods between the factory and the distribution center, which goal is to solve the location problem of the distribution center. The second level is to transport goods from distribution centers to customers, considering the capacity of the distribution center and customer demand and vehicle loading capacity, which is a classical VRP problem.
The establishment of different logistics distribution network systems needs to consider many factors, including whether to establish various positioning nodes, where to build and the number of positioning nodes. The cost includes start-up costs of facilities and distribution costs of vehicles. The objective function in this paper is to minimize location route problems, including the distribution center's start-up cost, transportation costs between distribution center to the customer, and vehicle start-up costs.
The VRP [19]-[21] is described by a directed weighted graph G: Given is the vertex set (0 is distribution center, others are customers). } , is weight matrix, ci,j is the distance between customer i and j.
There are most considerable capacity limits in vehicles, assuming all the vehicles are the same models. The research needs to arrange vehicles to service for all customers, and require meeting the following conditions and makes the vehicle shortest total route: • Any vehicle still cannot exceed its maximum load in the process of travel.
• Every customer is only serviced by one vehicle and served only once.
• Every vehicle starts from the distribution center and eventually returns to the distribution center.
is a solution of VRP, which is a vehicle route. VOLUME XX, 2017 1 The single factory multi-distribution center multicustomer center production model designed here is based on the following assumptions.
• Each customer is serviced only once, and the customer demand is uncertain.
• Each route starts and ends on the same distribution center, and the total load capacity of the route does not exceed the capacity of the vehicles service's route.

•
The total carrying capacity of each route should be no more than the capacity of the distribution center of serving the route.
Objective function (1) aims to minimize total cost including startup cost of distribution center, distribution route cost, and startup cost of vehicles. Constraint (2) requires each customer only belongs to one route, and only can be visited once. Constraint (3) ensures the possibility of vehicle capacity when demand changes. Constraint (4) ensures the possibility of the distribution center capacity when demand changes. Constraint (5) and (6) require each path begins and ends with the same distribution center. Constraint (7) ensures eliminating inner loop. Constraint (8) stands for if the vehicle through a certain path. Constraint (9) stands for if one distribution center is chosen.
The variable in constraint (3) and (4) is a triangular fuzzy number, which means fuzzy demand. For a given fuzzy demand , assumed that rest capacity of vehicle k is after delivering service nodes, and rest capacity of distribution center is by after servicing covered demand node. So and are also triangular fuzzy numbers: Then constraint (3) and (4) can be converted to the followings: This research uses a fuzzy c-means algorithm (FCM) [22]- [24] to generate the initial solution. FCM algorithm set c as parameters in the data set to divide n elements into c cluster. The group-dividing method is based on the membership degree of element and focus [25], [26]. The primary purpose is to make the weighted sum of squares is minimal in each data node to the cluster center. As shown in formula (10).
: Vector characteristics of customer node j. : Group center of clustering i. : Membership degree of customer node j to groups center of clustering i. m: Convergence parameter values, used to adjust the weight of membership degree, which is any number greater than 1, the number is higher, the convergence speed is faster. It is usually set to 2.
Membership degree can be calculated by harmonic mean. It's upper limit when = 1, which show customer node j is completely belongs to cluster i. Otherwise, it's lower limit when = 0, which show customer node j is completely not belongs to cluster i. = 0.5 shows the highest fuzzy degree. FCM has good clustering results and a low error rate and can deal with multi-dimensional data. Using the membership degree as clustering is more flexible and suitable for processing high volatility data, such as market segmentation, real-time customer market demand changes, medical cause analysis. Its main steps are as follows: Step 1: Initialization setting. Set minimum error values between the initialization group center and new group center and the update halt time.
Step 2: Decision group-dividing number. Group-dividing number must be decided before clustering. Set initial groupdividing number as the primary number of delivery vehicles, which can be obtained by dividing the total customer demand by each vehicle's car capacity.
Step 3: Set initial group center. Set c average as group center of the initial cluster after deciding group-dividing number c.
Step 4: Confirm membership degree. Firstly, calculates membership value of current customers. Secondly, check each customer's membership value. Thirdly, select t the highest membership value to assign the customer to the cluster. Judge whether according to the vehicle capacity constraints. By analogy, customers need to accumulate vehicle capacity to the cluster and move on to the next customer until all customers are members of the cluster.
Step 5: Determine whether there is any group center. Delete the vehicle of total capacity, and determine whether there is any group center. If so, then return to step 4. On the other hand, stop updating.

IV. ALGORITHM DESIGN
In the model proposed in this section, each supply chain member has its restrictive conditions in production, distribution, and inventory. When the supply chain members increase, its decision variables and restrictive conditions will increase substantially. If the exact algorithm is proposed to solve this problem, the solution process will not only be timeconsuming and difficult. In recent years, genetic algorithm and tabu search algorithm have been applied in this issue and acquired success. However, they have their own problems, such as tabu search algorithm has strong dependence for an initial solution, and the quality of feasible solution of genetic algorithm, on the whole, is not very high. The current research hot is a hybrid algorithm, by mixing to a certain extent, to make up for the defects of each algorithm.
The objective function can be divided into two parts: 1 ( ) f z is a location & transportation problem from factory to distribution center, which is a integer programming model, which can be solved by genetic algorithm.
2 ( ) f z represents each distribution center serves multiple customers, which is a typical vehicle routing problem. Tabu search algorithm is used to solve it. The GA-TS algorithm process is shown in Figure  1.

A. DESIGN OF GENETIC ALGORITHM
Genetic algorithm (GA) has been prevalent in solving optimization problems in recent decades because of their robustness in finding the optimal solution [27], [28]. The key to improvements is suitable mutation and crossover strategy [29], [30]. In this paper, an integer encoding method is adopted to generate pop size chromosome with a length of N (customer number) utilizing random generation, which is divided into vehicles according to the customer order of each chromosome according to the vehicle capacity constraint, and the initial population is obtained by decoding. Single genetic algorithm using serial number coding way, for a variety of products to multiple of the production scheduling problem, the production of the products by the order from top to bottom row from left to proper order, each representing a production plan. Because this paper assumes that only one production line, so it can be simplified as the serial number coding a 1 , a 2 , … , a n . The number of iterations to judge whether the set algebra, to choose the best performance after the chromosome of the corresponding solution as asked for production scheduling optimization solution of the problem output.

2) FITNESS FUNCTION
The objective function is the minimum cost of the distribution center location problem. The choice of the GA operations is to keep the individual of the largest fitness, so the objective function cannot be directly saved as an evaluation function. The reciprocal of the objective function is used as fitness function.

3) SELECTION STRATEGIES
There are many kinds of selection strategies, such as roulette, ranked selection, preserving the best individual, random league, graded selection, steady-state selection, etc. In order to ensure the overall nature of the algorithm, the sorting selection method was adopted in this paper to select the mother: all chromosomes in the population were sorted according to their fitness value from small to large, the sorting value of the chromosome with the smallest fitness value was 1, and so on, the sorting value of the chromosome with the largest fitness value was the total number of chromosomes. Then the ranking value directly replaces the fitness value as the probability in the roulette wheel method. The advantage of the sort selection method is that the selection is based on the relative position of the chromosomes in the population rather than the size of the fitness value so as to avoid the inability to highlight the better parent when the fitness value is too concentrated. The calculation steps are as follows: Step 1: Count the number of chromosomes N in the mother generation.
Step 2: Arrange all chromosomes in increasing order of fitness value.
Step 3: Give different sorting values t R according to their arrangement order, that is, give 1 to the chromosome with the minimum fitness value, and so on, give N to the chromosome with the maximum fitness value.
Step 4: Solve the probability of each selected chromosome.
Step 5: Solve the cumulative probability of each selected chromosome. TP t = ∑ p i t i=1 Step 6: Randomly generate N groups of random numbers S t between 0 and 1. If S t ≤ TP 1 , select chromosome 1; If TP 1 ≤ S t ≤ TP 2 , chromosome 2 was selected.
Step 7: Repeat step 6 until complete chromosome duplication.

4) CROSSOVER OPERATOR
In a regular crossover operation, the one-cut-node crossover operation is usually used by selecting a cut node randomly, where two parents originate two different offspring. In order to avoid possible generation of infeasible chromosomes after crossover operation, feasibility test should be carried out on offspring. On the premise of not destroying the population's genetic diversity, the algorithm's evolution speed is accelerated.
The crossover rate c P is adjusted adaptively according to the formula: Where, m f is the population's maximum fitness; ave f is average value of each generation's fitness. l f is the larger fitness value between the two individuals being crossed, 1 , 2 ∈ (0,1) . If 1 and 2 are set, can make adaptive adjustment.

5) MUTATION OPERATOR
The primary function of mutation in genetic algorithm is to increase population diversity and prevent entrapment into optimal regional solutions in the process of searching. The traditional mutation operation is to select one or more variables in the chromosome for mutation randomly. In this study, the modified multi-point mutation rule will be used as the basis for gene chromosome mutation. That is, several different locations will be randomly selected in the chromosome, and then the genes at this location will be calculated for mutation. The calculation steps of mutation rate m p are summarized as follows: Step 1: Randomly generate N groups of probability values t w between 0 and 1, where t = 1, 2, …N.
Step 2: t w compared with m p , if t m w p  indicates that the chromosome must undergo mutation operation.
Step 3: Determine whether the member needs to be mutated according to the value of 0 and 1 corresponding to the chromosome position.
If its value is 0, the member's gene value is retained. If its value is 1, it will determine which level in the supply chain it is located and then randomly generate purchasing decisions for each period according to its level.
Step 4: Repeat the above steps until the offspring chromosomes of group N are produced.
Step 5: Determine whether the mutant offspring chromosome conforms to each member's productive capacity or storage capacity constraints. If so, the chromosome shall be VOLUME XX, 2017 1 retained. Otherwise, the chromosome is deleted to avoid incomprehensible production.

6) STOP CONDITION
The stop condition of this study is that the evolutionary algebra executed has reached the set end algebra; then the calculation is finished.

B. TABU SEARCH ALGORITHM
Tabu Search (TS) [31], [32] is a meta-heuristic algorithm, which guides the local heuristic search process to explore beyond the local optimal solution space. Tabu search uses adaptive memory method to produce more flexible search behavior [33], [34]. In the calculation progress of the fitness function of GA, tabu search algorithm is used in a nested structure, which purpose is to calculate vehicle route cost in every optimization of selecting distribution. From the view node of raising the algorithm's speed, the research assign customers to their nearest distribution center and then use the following steps to solve VRP problems. This paper chooses the following algorithm strategy:

1) CONSTRUCTIVE AN INITIAL SOLUTION FOR THE TABU SEARCH
This encoding method directly produces an L arrangement of nonredundant natural numbers based on route sequence. Each customer is delimited to the distribution route according to the constraints of the distribution routing. The tabu rule is: if a customer node is chosen to insert, then it won't be selected again in the next generations until all the customers are included in the route, and so on. This meets the requirements of each node can only be visited once.

2) EVALUATION METHOD OF SOLUTION
For a distribution route of selected customer nodes, its advantages and disadvantages are determined by the objective function.

5) Termination criterion
Iterative operation end after a fixed number of algorithms, and the iteration length is set to 1000. The GA-TS algorithm in this paper firstly uses GA to conduct a global search. It makes individuals in the population more stably distributed in the solution space most areas and then starting from each individual in the group local search with TS algorithm to improve the quality of the group.
The better individuals generated by GA are arranged in increased order of cost, and the following TS procedure is performed for the top 60% individuals {Zc} in the list.
Step 2: Generate the initial solution Z' randomly from {Zc}, take the initial solution as the optimal solution Z' (Bsf = Z'), and calculate the fitness value F'(BestC = F').
Step 3: A new neighborhood solution set V {Z} = {Z1, Z2,.., ZTS} in TS is generated by randomly selecting the movement rules of the neighborhood solutions above. The fitness values are calculated to find the maximum fitness value F", and the corresponding solution is Z".
Step 5: From the neighborhood solution set, select the solution Z′′′ corresponding to the best fitness value F′′′ not in the tabu table as the optimal solution (Bsf = Z′′′, Best C = F′′′) and update the tabu list.

V. COMPUTER IMPLEMENTATION AND EXPERIMENTAL RESULTS
Numerical example has four distribution centers and a factory and 20 customers, which position are randomly selected from the coordinate's data in Solomon's test question bank of RC208. Node coordinate data are randomly generated within [0,100] 2 , and demand is randomly generated in [0,100]. The distance between each node uses Euclidean Distance of plane coordinates. ci,j = ¥0.8/km, 1 = 2 = 0.95. Each average startup cost of a distribution center is shown in Table 1. Each VOLUME XX, 2017 1 distribution center's start operation also needs to satisfy the maximum volumes, which are shown in Table 2. Table 3 shows each customer's position and fuzzy demand, Table 4 shows distribution centers' locations. Results are shown in Table 5 and Table 6.

All DCs 1 2 3 4
Code in GA 1 1 0 0  First, determine the appropriate population rules. For this purpose, the population size increased by 20 each time in the interval [60,300] test, 10 operations are performed for each population group. Figure 2 shows the average convergence of the algorithm algebra and convergence probabilities. As shown from Figure 2, with the increase of population size, the average convergence algebra generally decreases first and then increases, while the convergence probability is just the opposite. That is, too large or too small the population size will reduce the convergence probability. In particular, when the population size reaches 240, it has the maximum convergence probability and only needs the least number of iterations converges. The average convergence algebra is 48 and convergence rate is 0.09s. For this reason, the parameter of initial population size was selected as 240 is the most appropriate.

VI. DISCUSSION
According to the above parameters, this study conducted the experimental analysis with five sample data, respectively. The environmental characteristics of the sample data in the supply chain are shown in Table 7.
From the experimental results, the influence of the combination of genetic algorithm parameters on the solution is as follows: • When the parent's size was more significant, the quality of the solution was better, but the computation time was relatively longer. When the crossover rate was higher, the more chromosomes could perform crossover operation, the search space was larger, and the quality of the solution was relatively better.
• According to the experimental data of the crossover rate at three different levels, it can be found that when the crossover rate was 0.7, the solution obtained for most of the test examples was better than the solution obtained at the crossover rate of 0.5. However, compared with the solution obtained at the crossover rate of 0.9, it can be found that only the optimal solution of one example occurs at the crossover rate of 0.9. Therefore, a crossover rate of 0.7 was used for calculation in the experiment.
• When the mutation rate was larger, the situation was less likely to fall into the regional solution in the solution process. According to the experimental data of the mutation rate at three different levels, it can be found that when the mutation rate was set to 0.05, the quality of the solution was poor, and the variation of the solution was large. When the mutation rate was low, the solution had a higher chance of falling into the local solution, and there was no apparent convergence. When the mutation rate was set to 0.1, the quality of the solution was good, and the variation of the solution was slight. When the mutation rate was high, after 30 repeated tests, the solution values of each time converged more centrally. Therefore, a mutation rate of 0.1 was adopted for analysis. Based on the above results and considering the stability of the solution, the research adopted the following parameters: population size 240, and crossover rate 0.8, mutation rate 0.1. To verify the performance of proposed algorithm, the mathematical programming software LINGO was used to find the optimal solution. The genetic algorithm and tabu search were programmed by Visual C++. Figure 3 shows the results of 10 operations of the two algorithms. As shown from Figure 3, compared with LINGO, the iteration times of the GA-TS are about 30 to 90, which is much less than 23.7%. This indicates that the proposed algorithm is indeed effective. It has higher solved efficiency and can converge to the optimal solution faster. Table 7 shows the objective function values obtained by the two methods and the required solution time in. As can be seen in Table 7, using lingo and hybrid algorithms in this paper both can get the optimal solution. In terms of solving time, the use of a hybrid algorithm can significantly shorten the time, and the solving time of a hybrid algorithm will not be unstable with the expansion of the problem size.

VII. CONCLUSION
This paper stands for a strategic perspective of enterprise management and builds a production-distribution model in the supply chain. This is a highly complex NP-Hard problem, so the paper put forward the LRP model of the productiondistribution system and the algorithm design. A fuzzy caverage clustering algorithm is used to produce the initial solution. Fuzzy triangular number and confidence interval transformation are used to deal with fuzzy customer demand. Objective functions are decomposed and respectively solved by GA and TS. The experiment proves that this paper's optimization algorithm can effectively solve the LRP problem, which has a certain scale distribution center and fuzzy customer demand in the production-distribution system. The proposed hybrid algorithm of GA and TS algorithm overcomes the defect of local optimum in the literature viewpoint.
The research findings concludes that (i) determine the numbers of facilities with locations that should be opened in response to the quantity flow of products and (ii) minimize the total cost including startup cost of distribution center, distribution route cost, and startup cost of vehicles. This model can simulate single factory and single production, which also can simulate multi factory and multi-distribution centers and customers of the multi-level system in production-distribution plan. Future research can expand to real-time dynamic LRP problems, and more combination optimization methods can be tried to get a more optimum solution.