Performance evaluation of meta-heuristic algorithms for designing multi-objective multi-product multi-echelon supply chain network

The purpose of this paper is to develop and test three meta-heuristic algorithms to solve large size multi-objective supply chain network design problems. The algorithms are based on tabu search, genetic algorithm, and simulated annealing to find near-optimal global solutions. The three algorithms are designed, coded, and tested, and their parameters are fine tuned. The exact ε-constraint algorithm embedded in the General Algebraic Modeling System is used to validate the results of the three algorithms. The algorithms are compared using a typical multi-objective supply chain model utilizing several performance measures. The measures include the mean ideal distance, diversification metric, and percent of domination, inverted generational distance, and computation time. The results show that the tabu search algorithm outperformed the other two algorithms in terms of the percent of domination and computation time. On the other hand, the simulated annealing solutions are the best in terms of their diversity. The work in this paper is expected to help managers to solve large-scale supply chain problems that arise in oil and gas, petrochemical, and food supply chains.


Introduction
A supply chain network (SC) is a collection of information, workers, facilities, and processes necessary to turn raw material into finished products and transport them to the end user. Usually, all components of supply chain such as planning, production, marketing, and transportation have worked independently. Each part of the supply chain has its objectives, which are often contradictory. The integration of various goals into a supply chain system is a major issue.
Besides, real problems of these systems are complex and have non-commensurate and multiple goals in large nonhomogeneous decision space. Designing such supply chain networks with several echelons (stages), and with each echelon having several entities, is an important and challenging problem that has attracted and intrigued both practitioners and academicians. Most of the supply chain design challenges involve complex optimization problems. This type of problem is NP-hard. In past decades, these types of problems were formulated as single-objective optimization models to either maximize the profit or minimize the cost. Recently, in addition to economic objectives, organizations are interested in other objectives such as minimizing the negative impact on the environment, maximizing service levels, and minimizing risk. This shifts the supply chain network design problem formulation to a multi-objective optimization model and changes the optimality to Pareto optimality. For instance, Chavez et al. (2020) addressed a multi-objective supply chain problem integrating three objective functions, which are social impact, environmental impact, and operational cost. Kadzinshi et al. (2017) also proposed a multi-objective green supply chain system with three objective functions: carbon dioxide emissions, costs, and fine dust emissions. In addition, Hasani et al. (2021) integrated environmental and economic concerns into the design of a resilient supply chain system. In this paper, a typical multiobjective model is presented in Sect. 3 and is used to compare the performance of the developed algorithms.
The solution of multi-objective optimization models formulated for the design of multi-objectives supply chain network problems via exact algorithms such as the e-constraint algorithm is a time-consuming task and is impossible for large size problems. This necessitates the use of metaheuristics to find near global-optimal Pareto solutions for large problems in an efficient manner. The purpose of this paper is to develop three efficient metaheuristic algorithms based on tabu search algorithm (TSA), genetic algorithm (GA), and simulated annealing (SA) to solve large size multi-objective supply chain network design problems. The algorithms are designed and coded, and prior to testing them, their parameters are set properly to obtain robust and reliable solutions with high performance under various conditions. A regression model is developed for the algorithms' parameters setting. The results of the developed algorithms are validated using a single objective supply chain version of the model utilizing the results of the e-constraint algorithm embedded in the General Algebraic Modeling System (GAMS). This paper also provides a multi-objective and multiple products supply chain model to guide and help decision-makers in developing their strategic decisions within the framework of multiple dimensions. The provided model in this paper integrates a set of suppliers, plants, warehouses, distribution centers, and the end market/customer. Integrating all stages of the supply chain into one model leads to an integrated supply chain management that can be adapted to many practical situations that arise in real life as needed.
Well-designed experiments using a typical multi-objective optimization model for the design of supply chain networks described in the above paragraph are utilized to compare the performance of the developed meta-heuristics algorithms. The comparisons also included the exact econstraint algorithm. The algorithms were compared using several performance measures that include the percent of domination, diversity, and computation time. Sound statistical tests are implemented to differentiate between the algorithms. The results indicate that the tabu search outperformed all algorithms in terms of percent of domination and computation time while the simulated annealing solutions are the best in terms of their diversity. The nomenclature of the expressions used is shown in Table 1.
The rest of the paper is organized as follows: the literature review is provided in Sect. 2; in Sect. 3, the supply chain network design problem is described, together with a typical mathematical model for the problem. The design of the algorithms is outlined in Sect. 4, followed by data generation and setting the algorithms' parameters in Sect. 5. The comparison of the algorithms and their results is presented in Sect. 6. Section 7 concludes the paper and suggests directions for further research.

Literature review
Designing SC network problems has features of complex optimization problems. For finding the optimal solution, all possible combinations of the variables must be discovered. Dealing with the solution of these problems using exact algorithms consumes a great deal of time and may even be impossible for large networks. This motivates the use of meta-heuristic algorithms for designing supply chain network systems. The main advantage of meta-heuristic algorithms is the ability to provide acceptable near-global optimal solutions in a reasonable time. These algorithms examine parts of the search space to enhance the solutions and therefore reduce the time needed in finding a solution. Yu et al. (2022) conducted an extensive literature review to explore green technologies and sustainable supply chain management. Bibliometric networks were constructed  Khalilpourazari et al. 2020;Houssein et al. 2020;Kaveh & Zaerreza 2020). The literature review in this section focuses on metaheuristics algorithms for solving multi-objective supply chain (MOSC) design problems. One of the most common techniques is using genetic algorithms. The results obtained by using this techniques have proven that genetic algorithms can find near-optimal solutions. Altiparmak et al. (2006) utilized a GA to obtain efficient solutions for a mixed integer nonlinear programming (MINLP) model. In this model, long-term planning decisions of the supply chain activities including plant selection, collection center location, and distribution center sites and the quantities that should be shipped between entities that satisfy customer's requirements are determined. Lee et al. (2007) formulated a mixed integer linear programming (MILP) model for designing CLSC to minimize the total cost and maximize customer satisfaction. The model was transformed into a single objective using fuzzy goal programming approach. A GA was proposed to solve the model efficiently. He et al. (2007) proposed a novel hybrid method based on fuzzy approach and GA to optimize a deterministic reverse supply chain model. The design decisions are to establish the optimal location of solid waste and determine the number of treatment facilities. Xu et al. (2008) proposed a stochastic model tackling uncertainties in a customer demand and shipping cost. The authors formulated the problem as a stochastic MINLP model and then transformed it into deterministic MINLP. They considered strategic decisions such as selecting the entities, and the flow of quantities that meet the customer's demand. A novel approach based on a GA was used to design the network. Economic and customer service objectives are optimized. Farahani and Elahipanah (2008) presented a multi-periods, multi-products MOSC model. The objectives are to minimize the sum of surplus, backorders of products, and the total cost under supply capacity delivery and on holding constraints. A GA was used to obtain the solutions of the proposed model. Bevilacqua et al. (2012) developed an approach for distribution networks strategic configuration through nonlinear multi-objective optimization based on MOGAs. A real distribution chain of a large enterprise model is solved by the proposed genetic optimization.
Moreover, Sadigh et al. (2013) integrated a deterministic MOSC model for supplies selection and location of distribution centers. Quality, delivery time, and total cost are three objectives that are considered. Pareto optimal solutions are obtained using a modified GA. The results demonstrated the effectiveness of GA in obtaining solutions with high quality and reasonable time.
Mirakhorli (2014) considered and analyzed a model for a multi-objective closed loop supply chain (CLSC). The design model was considered under fuzzy environment to answer many strategic and tactical decisions such as locations of plants, collection centers, and distribution center sites and the quantities that should be shipped between these entities. A heuristic algorithm based on GA was employed in obtaining efficient solutions for this model. Pasandideh et al. (2015) incorporated uncertainty of demand, production and setup times, and costs into an MINLP supply chain model. The aim is to determine the production quantities, the flow through supply chain network, the shortages and inventories in each entity, and select the entities to minimize the expected and variance of the total cost. A GA was used to solve this model. Azadeh et al. (2017) integrated multi-objective midstream and upstream of an oil supply chain. A particle swarm and genetic algorithm were used to find efficient solutions of the model. Genetic algorithms are also employed for designing MOSC problems in Sahu (2013), Validi et al. (2014), Dzupire and Nkansah-gyekye (2014), Bandyopadhyay and Bhattacharya (2014), Arabzad et al. (2015), Yang et al. (2015). Furthermore, only a few studies used a tabu search algorithm (TSA) for supply chain problems. Shiguemoto and Armentano (2010) considered a multi-period, production-distribution model. They addressed two problems over a finite time horizon. Khalaf et al. (2011) formulated a MILP model for the design of the product family supply chain to minimize the total cost. The authors utilized a tabu search algorithm to minimize the total cost. A multi-echelon network problem was redesigned using a deterministic MILP model developed by Melo et al. (2012). To minimize the overall cost, various decisions are taken into account, such as the relocation of the entities, the volume of goods to be shipped between the entities, the investment budget that is available, and the inventory levels at storage facilities. A tabu search algorithm is used to solve this model. Mohammed and Duffuaa (2020) also developed a TSA for the design of MOSC systems. Previous studies have shown that the TSA is able to deal with combinatorial problems efficiently.
In addition, a few researchers employed simulated annealing for designing these types of problems. Mansouri (2006) utilized SA to solve a multi-objective sequence problem in a multi-stage supply chain. The model decisions are associated with organizing the set-ups stages of the supply chain. The results have shown the effectiveness of the algorithm for obtaining acceptable solutions.
On the other hand, a particle swarm is utilized for solving supply chain design problems for multi-objective cases. Shankar et al. (2013) addressed strategic and tactical decisions for designing a multi-echelon network. The decisions are to determine the production quantities, the flow through the network, and the number of facilities. The aim is to minimize the transportation costs and meet customer needs. A multi-objective hybrid particle swarm was utilized for finding the optimal solutions. Govindan et al. (2014) addressed an integrated sustainable multi-objective food supply chain model to optimize the number of facilities and the shipment vehicle size. Different algorithms including particle swarms were utilized to tackle this problem. A random-key approach was utilized for obtaining the solutions for continuous variables. Canales-bustos et al. (2017) addressed a multi-objective mining supply chain system. Economic aspects, along with environmental and emissions, are the considered objectives. A particle swarm was employed for solving the model. A particle swarm algorithm is also utilized in designing MOSC models (Mastrocinque et al. 2013;Park and Kyung 2014).
Moreover, a few researchers utilized an ant colony algorithm for the design of optimal supply chain problems. Sun et al. (2008) utilized an ant colony approach for dealing with a multi-objective supply chain design network. The findings indicated that the proposed technique is useful for solving such problems. Moncayo-Martínez and Zhang (2011) developed an efficient algorithm based on an ant colony algorithm (AC) to minimize the cost of supply chains and lead times. Other studies that utilized metaheuristic approaches for designing MOSC problems are conducted by Eskandarpour et al. (2013). The authors addressed a sustainable supply chain model to minimize tardiness and optimize cost and the environmental aspects. An efficient meta-heuristic algorithm based on a variable neighborhood search is proposed for solving the model. Devika et al. (2014) also presented several hybrid heuristics for addressing multi-objective closed-loop supply chains. These algorithms are based on an imperialist competitive algorithm and a variable neighborhood search algorithm. Sarrafha et al. (2015) developed a novel algorithm based on a biogeography optimization approach for the design of the MOSC system. The authors considered different decisions such as determining the production and flow of quantities and inventory of products at factories and distribution centers. A heuristic algorithm based on biogeography was used to solve this problem. Saffar et al. (2015) formulated a stochastic MOSC model for designing a green supply chain. Environmental and economic issues were considered as the objectives of this model. They applied a multi-objective differential evolutionary algorithm to solve the proposed model. Zhang et al. (2016) developed a bee colony algorithm for solving an MINLP supply chain design model. Kayvanfar et al. (2017) developed an algorithm based on water drops to solve a multi-objective supply distribution model for industrial clusters. Kadzinski et al. (2017) considered a MOSC model. The total cost along with fine dust and carbon dioxide emissions are the considered objectives. They utilized different evolutionary algorithms to tackle the developed model. Other work that used meta-heuristic approaches for designing MOSC problems were conducted by Pishvaee et al. (2010), Fallahtafti et al. (2014, Yu and Goh (2014), Govindan et al. (2015), Fahimnia et al. (2015), Nooraie and Mellat Parast (2015), Zhang et al. (2016).
In addition to the metaheuristic algorithms, other techniques are used to tackle the multi-objective issues. Alçada-Almeida et al. (2009) developed an interactive decision support system by combining different methods and goal programming to solve a multi-objective facility location model. Goal programming, fuzzy, and weighted sum multi-objective programming methods are used by Galante et al. (2010) to analyze the solution space. First, using goal programming, the objectives' value is assessed. Then, using fuzzy multi-objective programming and weighted sum techniques, a Pareto-optimal solution is found between these solutions. Pati et al. (2008), Ramudhin et al. (2010) also used goal programming. Ionescu and Vernic (2021) presented a new multi-objective algorithm based on symbiotic organisms search for maximizing the quality and minimizing the total duration and total cost of scheduling problems. Khan et al. (2021) evaluated green supply chain management strategies for the food supply chain. The authors used a hybrid method based on a fuzzy approach and distance-based assessment technique for comparing five strategies. In addition, Khan et al. (2022) used a simple sampling method to construct a conceptual framework based on the environment, organization, and technology theories to explore the relationship between a sustainable supply chain and a sustainable development strategy. Other techniques used in supply chain field include Lagrangian Decomposition Scheme (Wu and Golbasi 2004); variational inequalities (Yang et al. 2009); modified projection method (Nagurney 2009); ANP methodology and MIMOP method (Wu et al. 2009); L-shaped algorithm (Cardona-Valdés et al. 2011). Table 2 classifies the reviewed papers based on two main categories. The first one is based on the type of metaheuristic algorithm, whereas the second is based on the objective functions (cost, profit, risk, or others). Based on the classification; the majority of the published papers have employed genetic algorithms and particle swarms for solving the MOSC network design problem. In contrast, the literature indicates that there is a lack of efficient approaches for solving MOSC network design problems. In Overall, the literature review shows that the simulated annealing and tabu search are two meta-heuristic approaches that incorporate mechanisms for going beyond local optima. They have an ability of obtaining robust solutions for solving complex problems in many areas, such as the routing problem (Martínez-Puras and Pacheco 2016), the wireless networking problems (Huertas and Donoso 2018), problem of the allocation of resources (Belfares et al. 2007), the location of the facilities problems (Sun 2006;Li et al. 2009), sequencing problem (Mansouri 2006). Accordingly, it is valuable to develop simulated annealing and tabu search and investigate their performances for the design of MOSC networks. Moreover, the literature review indicates that the most used approach for solving the problem under consideration is the genetic approach, while The multi-objective scope is adopted for modeling the SC because of its versatility in providing various trade-offs between possible solutions Fig. 1 Proposed supply chain network the simulated annealing and tabu search are not utilized to the best of their potential. In addition, no study has compared the performance of the three approaches (GA, TSA, and SA) in solving supply chain network design problems integrating profit, risk, and environmental aspects in a single model considering all stages of the supply chain network. This paper will address this research gap by developing three algorithms: one based on TS, the second one based on GA, and the third one based on SA. Then, a comparison was made of their performance for solving different sizes of the supply chain network design problem using several different performance measures and sound statistical tests. In short, the advantages, disadvantages, limitations of the proposed algorithms, and novelty of the mathematical model of the current paper are shown in Table 3.

Problem description and the proposed mathematical model
The large SC network depicted in Fig. 1 is taken into consideration in this section, which consists of suppliers, factories, warehouses, distribution centers, and end users. The raw materials are supplied by suppliers, shipped to factories, and then converted into useful products in the factories. After that, the products are shipped to warehouse centers that transport the products to the distribution centers (DCs). The distribution centers store the products before they are delivered to the end customers. The problem is to design a supply chain network that allows for strategic decisions such as optimizing the number and location of suppliers, plants, warehouses, and distribution centers to be established, and operational decisions such as optimizing the flows through the network. Three objectives lie in the considered multi-objective model. The first goal is to maximize the total profit of the supply chain in order to give shareholders as much wealth as possible. The second objective is to minimize the supply chain risks, although natural crises are the main sources of risk, other risk sources are associated with the poor quality of raw materials and products and the delay in the delivery time of items. Therefore, this work tries to minimize the quality and delivery risks of the supply chains. The third objective is to minimize the supply chain emissions. Because of the growth in awareness regarding climate change, companies are enforced to reduce the emissions from doing business. The proposed model aims to achieve environmental sustainability. The assumptions of this model are: 1. The integrated supply chain is designed within a fixed planning horizon.
2. Potential locations and the capacity of all entities are known in advance. 3. Customers' demands are predetermined and known. 4. The capacities of all the warehouse and distribution centers are fixed and limited for storing products. 5. The production capacities of the manufacturing plants are fixed and limited.
The formulation of the supply chain problem is described above and depicted in Fig. 1

Sets
S Set of indices (s = 1,2, … ,S) used for suppliers I Set of indices (i = 1,2, … , I) used for potential factories J Set of indices (j = 1,2, … ,J) used for potential warehouses K Set of indices (k = 1,2, … ,K) used for potential distribution centers M Set of indices (m = 1,2, … ,M) used for demand markets P Set of indices (p = 1,2, … ,P) used for products T set of indices (t = 1,2, … ,T) used for raw materials Parameters EP i Unit environmental effect of producing products at factory i;

ETS sit
Per unit distance environmental impact of shipping raw material t between supplier s and factory i;

ETP ijp
Per unit distance environmental impact of shipping product p between factory i and warehouse j; ETD jkp Per unit distance environmental impact of shipping product p between warehouse j and DC k; ETW kmp Per unit distance environmental impact of shipping product p between DC k and market m; Prd st Delivery risk probability for raw material t supplied by supplier s; The model optimizes three objective functions. The first one, Z 1 is to maximize the total profit, which is represented by Eq. (1). The total profit is the difference between total revenue and total supply chain cost. The total supply chain cost is a combination of the fixed cost of establishing the plant, warehouse, and DCs, the purchase cost of raw materials, the production cost at plants, and transportation cost of shipping quantities through the supply chain network. The second objective, Z 2 , is to minimize the expected supply chain risk formulated as shown in Eq. (2). This paper defines the expected supply chain risk as the expected quality and delivery risks caused by supplying the raw material from suppliers, the expected quality caused by producing the product at plants, and the delivery risks caused by shipping the products through the supply chain network. The third objective, Z 3 is to minimize the supply chain emissions, which is represented by Eq. (3). The emissions are caused by establishing plants, warehouses, and distribution centers facilities as well as the emissions caused by producing products at plants, and the emissions due to transportation. Max Min In this paper, the model optimizes the three objectives described in the above paragraph while satisfying many practical constraints represented in Eq. (4) through Eq. (13). Equation (4) represents customer's demands constraints. The customers' demands for each type of product must be satisfied. Other constraints are material balance which are represented by Eqs. (5-7). The third type of constraints is the capacity constraints of processing entities. These types of constraints are represented by Eqs. (8)(9)(10)(11). Each facility at each echelon of the supply chain has a limited capacity. The total quantities of the products transported to a facility must not exceed its capacity. The nature of decision variables is represented by QSI sit ; QIJ ijp ; QJK jkp ; QKM kmp ! 8s 2 S; i 2 I; j 2 J; k 2 K; m 2 M; p 2 P; t 2 T 4 Proposed algorithms Many combinatorial techniques are utilized for solving large-scale optimization problems. These approaches have an ability to search the space of the problem intensely and then solve the problems efficiently and effectively. Accordingly, three efficient meta-heuristic algorithms are proposed to find near global optimal solutions. The developed algorithms are based on tabu search, genetic algorithms, and simulated annealing. To do this, the required structures of the meta-heuristic algorithms are developed. In this paper, although a set of constraints including continuous variables define the search space of the proposed problem, the binary variables define an attractive feasible space. This property is exploited to develop an efficient solution approach for solving such problems. In this paper, the problem is divided into two parts. The first part is to set the optimal strategic SC network design by establishing the best location for the facilities (plants, warehouses and distribution centers) in each echelon and selecting the best suppliers. Once the structure of a supply chain network is set, the second part is to determine the short-term decisions such as the production quantities and the shipments of products through the supply chain.
In the first part, the supply chain design structure is selected using the proposed meta-heuristic algorithms. Once the facilities are selected, the second part of the problem is handled by a linear programming method. For each structure of the supply chain design solution acquired in the first part, the optimal flow through the supply chain network is obtained in the second phase. Reducing the computational time and the search space are the main benefits of this approach. In addition, the proposed approach can be utilized to solve the high-dimensional problems intensely through searching the problem space comprehensively.
The facility is chosen to be open at random as part of the diversification procedure in this paper. Instead of emphasizing the solutions' quality, this process emphasizes their variation. The solutions for the supply chain network design structure that cannot satisfy the required demand are removed for reducing the computation time. In this regard, a condition to check the feasibility of open facilities in each echelon is proposed. The open facilities must be able to satisfy the demands. In this paper, a random weight method is utilized to find the most efficient solutions (Miettinen 2012;Murata et al. 1996). This approach gives all possible efficient solutions a uniform chance and helps to explore the whole search space. The m objectives are converted into a single objective function by assigning a random weight w i for each objective as follows: . . .; m and X m i¼1 w i ¼ 1 In this paper, w i is calculated as follows: where N i is a random positive number, N i 2 1; . . .; 100 f g : Due to the difference in the units of the objectives, the formula in Eq. (16) is applied to convert them to non-dimensional quantities (Pishvaee et al. 2010).
In case objective i needs to be maximized, the inverse ratio is used.

Multi-objective tabu search (MOTS) algorithm
The proposed algorithm used in this paper is based on tabu search algorithm developed by Glover (1990). The developed algorithm has the ability to acquire optimal and nearoptimal solutions within reasonable amount of time. In the proposed algorithm, the purpose is to obtain good solutions with a sufficient diversity to cover the entire Pareto frontier within an acceptable time. Various exploitation and exploration strategies are employed. The tabu list is the first element of the proposed tabu search. The solutions that have been visited are saved in a tabu list; thereby, revisiting these points is prevented. The allowable number of solutions kept in the list is fixed. In addition, an aspiration condition is used to rescind a tabu transition in order to improve the current solution. Another principle is the structure of the neighborhood. This principle aids in the transition from one solution to another. In this paper, the type of neighborhood structure is based on changing the status of a facility. In order to reduce the computation time for obtaining an optimal solution, only a fraction of neighbors in the neighborhood are evaluated. This could also act as a diversification procedure. During the algorithm search, a backtracking mechanism is applied to avoid local optima. The intensification procedure is another tactic used in the development of a tabu search algorithm structure to enhance the quality of the achieved solutions. The tabu search algorithm is restarted with the best solution after a limited number of iterations. The backtracking methodology and the intensification phase are only used once during the process of algorithm search. In this paper, the initial heuristic algorithm, strategy of local search, and other proposed approaches are combined with the tabu search algorithm to construct the developed solution approach. The pseudo-code of the tabu search algorithm is shown in algorithm 1.

Multi-objective genetic algorithm (MOGA)
In this paper, MOGA is developed to solve the MOSC design problem and the results obtained are compared with the other proposed algorithm results. In this paper, the MOGA developed by Altiparmak et al. (2006) is used. The first step for the successful implementation of a genetic algorithm is designing the chromosome properly because the chromosome representation plays an important role that effects the quality of the solutions. In this paper, the strategic design variables represent an appealing search space in the considered problem. This gives an opportunity to represent the chromosome as a single dimensional array. The array includes binary values that represent strategic decision variables related to the suppliers, plants, warehouses, and distribution centers, respectively, as shown in Fig. 2. Algorithm 2 shows the pseudo-code of MOGA.
In this section, the initial population is generated as follows: first, an initial solution is generated using algorithm 4. Then, 2-exchange, 4-exchange, and swap operators are applied to construct an initial population. The aim of using more than one neighborhood structure is to control the search intensification and diversification. The fitness function value for each solution (chromosome) is calculated using Eq. (14).

Crossover and mutation operators
In the GA, the children are generated by performing an operator called a crossover. This operator helps to discover a new search space through mating two chromosomes. A segment-based crossover operator is employed. In our problem, each chromosome has four segments that represent suppliers, plants, warehouses, and distribution centers, respectively. For each segment, binary numbers are assigned with equal probability for all the segments. This procedure is called a binary mask. Figure 3 shows the crossover and the binary mask, while ''0'' means that the gene of the child is transferred from the first parent and ''1'' means that the gene of the child is transferred from the second parent. The characteristics of both parents are transferred to the new child through the crossover. The crossover operator is applied probabilistically to a pair of chromosomes in the current generation to generate offspring for the new generation with crossover probability (Pc). To reduce the computation time, one of the new children is selected to be survived based on the fitness function. The roulette wheel mechanism is utilized to select a pair of chromosomes.
On the other hand, the mutation operator helps to improve the diversity of the new generation. It is performed by changing the chromosome gene. Similar to the crossover, a segment-based mutation mechanism is used as shown in Fig. 4. In this operator, a binary random is first assigned randomly for each segment of the chromosome. Then the selected segment is mutated. In the mutation operation, a two-exchange operator is utilized by changing the status of the two genes from the selected segment. The two genes within the segment are selected randomly.

Selection mechanism
In the proposed genetic algorithm, the initial population is constructed as described in the above section and the nondominated solutions set is constructed according to the dominance optimality condition. At each generation of a genetic algorithm, the non-dominated solutions set is updated. In the selection mechanism, the survive chromosomes for the next generation are selected based on the (npop ? q) selection strategy, where npop is the population size and q is the size of the children in the current generation. First, two solutions from the efficient solution frontier are selected randomly and are then used to construct the new generation. After that, (npop -2) different best chromosomes from the constructing evolving pool are selected. The weighted sum method is utilized to identify the best chromosomes.

Validation
In the proposed GA, only the binary variable solutions that satisfy the customer demands are evaluated. In this regard, a condition is added to validate the solution's feasibility by checking the capacity of the binary chromosome. If the Performance evaluation of meta-heuristic algorithms for designing multi-objective multi-product… 12233 feasibility condition is not satisfied, crossover and mutation operators will be backtracked and the check performed again. The idea behind the feasibility condition is that the total capacities of the open facilities at each echelon must satisfy the total customer's demands.

Multi-objective simulated annealing (MOSA) algorithm
In this paper, a simulated annealing algorithm is also proposed for solving the MOSC design model. Kirkpatrick et al. (1983) developed the first version of simulated annealing for solving complex problems. In the simulated annealing algorithm, the temperature and cooling rate play an important role in determining the properties of the material structure. First, the system starts with a high temperature, and then, the temperature is decreased gradually until the system reaches a stable state. The idea behind this is to avoid local optima and search for global solutions. In developing a simulated annealing structure, the initial feasible solution is generated using algorithm 4. The structure of the solution is similar to the structure of the MOTS. The search moves from a certain solution to another one by changing the status of a facility. If the generated solution dominates the current solution, the current is replaced by a new solution. Otherwise, a given probability is assigned to accept the generated solution. The probability of acceptance (P(A)) is provided by the following formula: where DE is the difference between the weighted objective functions of the current and new solutions, and T is the temperature. The temperature is constant during a certain iteration number (K). After that, the temperature is decreased by a cooling rate (a), where a 2 ð0; 1Þ Abido (2000). Algorithm 3 shows the pseudo-code of the developed simulated annealing. In this paper, the initial temperature levels are also selected based on several preliminary experiments to accept the solutions with a probability of more than 0.95 at the first iterations and to accept only the best solutions at the final iterations.

Construction of an initial solution
A good initial solution usually reduces the running time of meta-heuristic algorithms and influences the quality of the solutions. A competitive heuristic algorithm is proposed in this regard to constructing a feasible starting solution. Two stages are used in developing the heuristic algorithm. The first one is to set the structure of the supply chain network design by selecting and establishing the facilities among the available locations at each echelon. On the other hand, the flows of the products through the supply chain network structure that are set in the first stage are optimized in the second phase. The binary variables are selected based on the lower combination of the fixed and the transportation costs using the proposed formula (18).
where f 0 i is the financial loss at the suppliers and represents the fixed cost of establishing facility i in the rest of the supply chain echelons, C ij represents the shipping cost of the flows between facilities i and j, and b i represents the total capacity of facility i. Based on Eq. 17, the facilities that have a lower cost are open as necessary to satisfy the demand. Once the binary variables are selected, the flows through the network are optimized using a linear programming method. Algorithm 4 shows the pseudo-code for getting a starting solution that is both acceptable and feasible.

Performance evaluation
To evaluate the performance of the proposed metaheuristic algorithms, several computational experiments are designed. To do this, the performance metrics of multiobjective problems are presented. Then, the parameters of the proposed algorithms are set using a regression analysis technique. Finally, the results of the proposed algorithms and e-constraint method embedded in GAMS are compared. Several statistical tests are used to conduct the comparisons. However, the e-constraint method embedded in GAMS obtained solutions for small and medium size instances. The pseudo-codes of the proposed algorithms are developed and built using MATLAB R2017a. The problems are solved on 2.6 GHz computers with 4 GB RAM.

Evaluation metrics
The quality of efficient solutions obtained by the proposed algorithms is measured by three important characteristics, which are diversity, convergence, and computation time. Consequently, in this paper the most common metrics that measure the three features that are used to fine-tune the parameters of the developed algorithms are: 1. Mean ideal distance (MID) (Behnamian et al. 2009): evaluates the closeness of solutions with respect to the ideal solution. It is given as follows: where c i ¼ F i À F ideal k k , n is number of Pareto solutions, and F ideal is the best solution in the case of solving the problem as a single objective. 2. Spread of non-dominance solution (SNS) (Behnamian et al. 2009): to evaluate variation in solutions from the best points. It is calculated as follows: The algorithm with higher SNS value is preferred.
3. Diversification metric (DM) (Deb and Jain 2002): to measure the solution range. DM is calculated by: where max i f ji is the maximum value of the objective function j obtained by the algorithm i, and min i f ji is the lowest value of objective function j acquired through algorithm i.
The algorithm with higher DM value is better.
4. Percent of Domination (POD) (Zitzler and Thiele 1998) evaluates the dominance rate of the solutions obtained by algorithms. To compute this metric, first, all the solutions obtained by the algorithms are mixed. Then, the dominated solutions are eliminated and the ratio of the efficient solutions belonging to each algorithm is computed. It is computed as follows: where Z 0 Y means that the solution Z dominates solution Y, alg i is the algorithm i, and P i is the Pareto solutions number of algorithm i. The higher POD value is better. 5. Inverted Generational Distance (IGD): this metric is developed by (Van Veldhuizen andLamont 1998, 2000). It refers to the distance between the acquired Pareto solutions and the true Pareto front (PF true Þ. IGD can be calculated as follows: where A is the set of near-optimal solutions obtained with the algorithm, z is the Pareto front element, a is a point in A, and d is the Euclidean distance. 6. Computation time

Data generation
To assess the performance of the developed algorithms and to set their parameters, several random problems are generated and tested. The test problems are divided into three levels as on the facilities number at each echelon, namely small size, medium size, and large size as shown in Table 4. The parameters for the supply chain networks are provided in Table 5. Other parameters of the supply chain network design are computed based on the given parameters to keep them in a realistic structure.

Setting algorithms parameters
Proper tuning of the algorithm parameters leads to enhancing the solution quality. Therefore, the developed algorithm parameter values are tuned properly. First, different algorithm's parameters levels are chosen based on the size of the problem, conducting several primary analyses, and the literature as shown in Table 6. The algorithms on the medium problem of each scale are employed 36, 27, and 27 times for TSA, GA, and SA, respectively. For instance, three parameters need to be set for the tabu search. These parameters are tabu size (4 levels), neighbors percentage (3 levels), and number of iterations (3 levels), resulting in a total number of combinations of 36 times. The medium problem is selected to reduce the number of experiments and run times. Furthermore, it is impractical to optimize the algorithm parameters for each test problem because it requires a long time. Although the algorithm parameters are set based on a medium scale, some optimal parameter values are determined by the size of the problem. Each time, the values of the algorithm parameters are changed within corresponding levels to acquire the performance metrics. Moreover, under a fixed combination of the algorithm's parameters, each algorithm is employed for ten runs under different random environments (e.g., weight w i , x). Then, the averages of the runs for the medium problem of each scale are computed.
In this paper, all of the algorithms' parameters are set using effective methods. In this regard, comprehensive statistical analyses are conducted to tune them properly. First, the relationships between the algorithm parameters and the performance metrics are described using a regression method. In this paper, the input variables are controllable by experiments and continuous with neglecting errors. A significant relationship between the inputs and response variable is estimated using MATLAB software. This relationship is represented by the quadratic regression function (Neter et al. 1996;Pasandideh et al. 2015) as shown in Eq. (24): where E Y ð Þ is the expected value of the performance metric, b 0 , b i , i = 1, 2, 3,b ii , i = 1, 2, 3, and b ij , i & j = 1, 2, 3, 4, i = j are the intercept, the linear term coefficients, the quadratic coefficients, and the interaction term  4 (8, 8, 8, 8, 20, 4) P.5 (10, 10, 10, 10, 40, 4) Medium scale P.6 (12, 12, 12, 12, 60, 4) P.7 (15,15,15,15,70,4) P.8 (18,18,18,18,80,4) P.9 (20,20,20,20,90,4) Large scale P.10 (22,22,22,22,100,6) P.11 (25,25,25,25,100,6) P.12 (28,28,28,28,120,6) P. 13 (30,30,30,30,130,6) coefficients, respectively, and X i ; i ¼ 1; 2; 3 are the algorithm's parameters. Afterward, the desirability function is utilized to optimize the values of the algorithm parameters within their corresponding levels. This method was developed by Derringer and Suich (1980) to optimize multiple responses. If the response has to be maximized, the desirability (d i ðb y i Þ) is described by the following equation: The equation for d i ðŷ i Þ, when it has to be minimized, is where s and t are the importance of each performance metric, b y i ðxÞ is the estimated value for response variable y, L i and U i represent the lower and upper bounds, respectively. In this paper, all metrics are equally important except the computation time in the large size problems. Subsequently, the values of the responses are converted into a desirable value that is calculated through the following equation: where n is the number of metrics. As a result, all values for the algorithm parameters, the R-Square, and the desirability for each proposed algorithm are shown in Table 7.

Results and discussion
This section discusses the results of comparing the proposed algorithms for solving multi-object supply chain problems. The comparison is based on several performance measures. In addition, several statistical tests are conducted to determine the best algorithm.  Performance evaluation of meta-heuristic algorithms for designing multi-objective multi-product… 12237

Comparison among the algorithms
The performance of the proposed algorithms for the design of MOSC is compared using several performance measures. The measures include percent of domination, solution diversity, and computation time. In addition, the proposed algorithms results are validated by the results of the e-constraint method. In this paper, we limited the number of solutions searched by the genetic algorithm and simulated annealing to the number of solutions searched by the tabu search algorithm. This procedure allows for comparisons of the algorithms using the same yard stick.
First, all problems are solved as a single objective using the developed algorithms and e-constraint method embedded in GAMS. The idea behind this is to validate the optimal values of the algorithm parameters. The values of each objective resulted in solving the problems as a single objective using the proposed algorithms are shown in Table 8. The results reveal that there is a slight difference between the results of the developed algorithms and the results of the exact method instilled in GAMS. This is an indication of the effectiveness of the developed algorithms with optimized parameters. It can also be observed that the tabu search algorithm obtains the best results among the proposed meta-heuristic algorithms when solving the problems as a single objective. This is due to the design and the characteristics of the developed algorithm, which guides the algorithm toward the path of global solutions and to search for more spaces, resulting in better solutions. After tuning the algorithms' parameters, the proposed algorithms are utilized for the design of the multi objective supply chain problems. Figure 5 illustrates the Pareto optimal solutions of a single run obtained by the developed algorithms at all scales. The result is also compared to the optimal Pareto solutions obtained by the exact algorithm embedded in GAMS as shown in Fig. 5a, b. The sample of solutions shown in Fig. 5 reflects the POD among the algorithms. As can be seen, the non-dominated solutions   Performance evaluation of meta-heuristic algorithms for designing multi-objective multi-product… 12239 Prob.8 Prob.9 Prob.10 Prob.11 Prob.12 Prob.13 MID MOGA 11,779,343 82,277,325 40,965,454 57,433,167 68,140,599 492,307,596.9 107,483,122.5 1,750,115,527 334,011,397 8,653,622,296 8,635,511,098 2,295,349,550 3,756,105,740 MOSA 11,796,552 81,250,448 39,126,236 57,670,326 69,676,241 491,994,175.6 106,951,415 1,745,654,966 334,061,261.4 8,650,395,541 8,632,347,219 2,309,156,624 3,763,516,129 MOTS 11,827,320 82,841,381 40,351,160 58,658,155 67,323,108 495,354,634.8 109,568,762 1,747,304,452 339,826,148.6 8,660,517,212 8,633,727,169 2,280,224,199 3,742,371,628 SNS MOGA 2,949,080 5,255,060 4,064,258 5,553,359 5,720,683 5,089,852.103 5,499,420.267 7,074,498.778 6,968,421.097 7,748,692.678 9,678,728.32 16,660,278.9 16,024,478.4 MOSA 3,222,456 6,583,532 4,711,693 6,175,994 7,300,760 6,362,972.645 6,110,407 9,399,708.931 8,442,487.439 12,708,951.76 14,684,047.8 18,104,638.1 16,686,729.4 MOTS 3,066,016 5,199,882 3,887,345 6,436,757 5,507,736 5,365,696.971 4,785,269 6,303,348.954 6,489,076.606 8,402,942.495 9,658,809.69 12,309,790.2 14,659,054 DM MOGA 9,033,200 18,269,420 16,853,480 21,386,630 22,951,222 21,803,230 23,813,817.5 29,346,237.42 30,102,789.72 32,815,000.1 39,779,476.4 58,248,452.16 62,870,446.38 MOSA 11,092,134 24,564,870 20,028,153 23,819,462 28,552,905 25,202,260 28,419,342 38,523,103.04 35,531,879.62 51,179,330.9 56,880,916 72,789,131.14 68,534,816.48 MOTS 11,082,464 21,560,936 18,230,510 27,755,389 24,196,302 24,966,663 19,817,256 29,305,200.19 31,064,301 36,548,517.3 40,431,311.1 56,424,676.55 67,209,162 obtained by the tabu search are closer to the GAMS solutions than the other two meta-heuristic algorithms. It is also noted that, as the problem size increases, the effectiveness of the MOTS is clearly observed. The results of problems 2, 9, and 10 are reported in this study to demonstrate the effectiveness of the proposed algorithm in term of POD. Furthermore, problem 10 is not provided when using GAMS because GAMS cannot solve such large problems. The GAMS computational time for solving this problem is more than 72 h without output. Afterward, the performance of the developed metaheuristic algorithms for the design of multi-objective, multi-product, multi-stage supply chain networks is assessed and compared based on five performance metrics (MID, DM, SNS, Computation time, and POD) as shown in Table 9. The performance measures are first computed based on the non-dominated solutions acquired by each algorithm. Since the MID, SNS, DM, and computation time values increase dramatically with increases in the problem Performance evaluation of meta-heuristic algorithms for designing multi-objective multi-product… 12241 size, the metrics values are normalized using the following formula: where f i algo ð Þ is the metric value of the algorithm i after normalizing, f i is the current metric value, and f min is the minimum metric value obtained by the algorithms. Table 10 shows the average of the metrics after normalization.
From Table 9, it is observed that the developed MOTS is superior to MOGA and MOSA in terms of POD and IGD in most of the instances. The average dominance ratio of the Pareto-optimal solutions of MOTS is 74.4% versus 39.9% for MOSA and 54.5% for MOGA, as shown in Table 10. On the other hand, the average IGD values for MOTS, MOGA, and MOSA are 1.933, 1.602, and 1.347, respectively. The key explanation for the excellent performance of MOTS is that the TSA features are useful to guide the algorithm toward a global solution direction and explore Moreover, as mentioned earlier, the intensification and diversification phases have a significant effect on the quality of the solutions and exploring more parts of the solution space. The values highlighted in bold indicate the best solutions. Moreover, it can be observed from Table 9 that simulated annealing has a lower ratio of domination than the other algorithms. The reason behind this is that the process for searching does not enter the most promising search space due to a lack of understanding the lack of the entire search space condition. With respect to computation time, the tabu search algorithm is superior to the MOGA and the MOSA. This is due to the fact that the tabu list prohibits many moves and prevents cycling. The computation times for the genetic algorithm are higher than the tabu search for all problems. The reason behind this is that the genetic algorithm has more operations such as crossover and selection.
Furthermore, the results show that MOSA is superior to MOTS and MOGA in the terms of the diversity of solutions. At the initial iterations of the simulated annealing, the algorithm allows for the acceptance of any moves even negative ones. The acceptance of such moves will become increasingly rare as the temperature falls. By accepting worse moves, the diversity of the algorithm will be improved. This results in higher diversity for the simulated algorithm than other proposed algorithms. Also, it can be observed that the convergence rate of the Pareto frontier members to the ideal point is 1.00979 for MOTS venues, 1.007 for MOGA, and 1.004541 for MOSA, respectively, as shown in Table 10. The reason that the simulated Performance evaluation of meta-heuristic algorithms for designing multi-objective multi-product… annealing has the best convergence rate is due to its ability to conduct good local searches. Moreover, numerous statistical tests for algorithm performance measures are performed in order to determine the best algorithm. Initially, the analysis of single factor variance is used to determine whether there is a statistically significant difference in all performance measures between the developed algorithms. Table 11 shows the results of a single factor of variance analysis for all metrics. Table 11 shows that the p-values for comparing the solutions achieved by the developed algorithms in terms of SNS, DM, computation time, and POD are less than 0.05. This indicates that the proposed algorithm solutions differ significantly.
Furthermore, additional statistical analyses, such as multiple comparisons of Tukey's assessments, are performed to determine the differences in detail (Montgomery 2009). MINITAB software is used in this paper to run Tukey's tests. Tukey's tests for all metrics taken into consideration are displayed in Table 12 after multiple comparisons. There are significant differences between the performances of the proposed methods, as indicated by the P-values less than 0.05. From Table 12, it is evident that there is no significant differences among the proposed algorithm solutions in terms of MID and IGD, whereas there are significant differences among the developed algorithms solutions in other terms.
Moreover, to support the results of the above analyses to compare the performance of the developed algorithms, relative percentage deviation index (RPD) is used (Khalilpourazari et al. 2016). As a result of this analysis, there is no significant statistical difference between the means of the measures in the case when the confidence intervals for the measure overlap. Accordingly, for each problem, the results acquired by the algorithms are converted into a relative percentage deviation index. This index is calculated as follows: where ALg sol is the metric value obtained by a given algorithm for a problem, and Min sol is the best solution obtained for that problem. The lower value of RPD is preferred. The means and the confidence intervals for the three algorithms are shown in Fig. 6. Referring to Fig. 6a, f, the MOGA, MOSA, and MOTS confidence intervals overlap for the MID and IGD measures, and hence, there is no significant difference among them with respect to these measures. Figure 6b shows the SNS measure, SA has the best performance among the three algorithms, and there is no significant difference between the other two algorithms. Figure 6c also demonstrates that SA is more effective than the other two algorithms with respect to DM. With respect to the computation time and the POD as shown in Fig. 6d, e, respectively, there is no significant difference between MOGA and MOSA for both measures. However, MOTS is significantly better than the other two algorithms for both measures. In conclusion, the MOSA has the most diverse solutions while the MOTS has the best performance for POD and computation time.
In short, Table 13 shows the proposed algorithms ranking in terms of the performance measures based on the Tukey's test results. The results in Table 13 show that the MOSA has a better performance comparing to MOGA and MOTS in terms of MID, SNS, and DM, whereas MOTS has the best performance among the three algorithms in terms of CPU-time, POD, and IGD.
The computation time for the metaheuristics and the exact e-constraint algorithm embedded in the GAMS are shown in Fig. 7. It can be seen that: the computation time of the developed metaheuristic algorithms is very low relative to the e-constraint method. The e-constraint method consumes a lot of computational time, which varies from 1.45 to 23.5 times the proposed met-heuristic algorithms in small and medium size problems. This result demonstrates the ability of the proposed metaheuristic algorithms for solving the MOSC design systems. It is worth noting that the exact algorithm was not able to provide solutions for large size problems. For example, the computation time of problem 10 is greater than 48 h. The tabu search algorithm has the smallest computation time among the three metaheuristics.

Discussion
This study proposed three metaheuristic algorithms for solving one of the common supply chain systems, multiobjective multi-echelon multi-product supply chain problem. The algorithms are based on tabu search algorithm, genetic algorithm, and simulated annealing algorithms. This paper contributed to the body of knowledge in the field of supply chain in several ways. At the model level, our model integrates all supply chain stages including suppliers, plant, warehouses, distribution centers, and end users. In addition, it integrates three objective functions: total profit, risk, and emission impacts. The findings of this study provide the decision-makers a list of optimal Pareto solutions. The decision-makers can select the best plan based on their preferences. Furthermore, the proposed solution algorithms can optimize any complex supply chain system and obtain near-optimal solutions within reasonable computation time. Since the quality of the solutions acquired by the algorithms depends on the initial solution, our study proposed an efficient heuristic algorithm to construct an initial feasible solution. The proposed algorithms' parameters are determined using a statistical approach based on regression analysis. This will enhance the quality of solutions obtained by the algorithms. The managers could use this methodology to solve large supply chain problems, analyze different solutions, and make decisions based on near-optimal values. The algorithms could be also applied in the future to any kind of similar models.
7 Conclusion, future research, and implications

Conclusion and future research
This paper develops three metaheuristic algorithms for solving the MOSC network design problem. The algorithms are based on the tabu search, genetic, and simulated annealing approaches. A competitive heuristic algorithm is developed to construct an efficient starting feasible solution for each of the algorithms. The algorithms are coded and tested and their parameters are fine-tuned using a quadratic regression approach. Then, their performance in solving the multi-objective supply chain network design problem is comparted using several performance measures. Sound statistical tests are implemented to evaluate the efficiency of the algorithms. The results demonstrate the effectiveness of the proposed algorithms for solving MOSC problems. Furthermore, the computational results indicate the superiority of the developed MOTS in comparison with MOSA and MOGA based on the percent of domination and computation time. On the other hand, the developed multiobjective simulated annealing has the most diverse solutions.
In this work, it is assumed that the production, warehouse, and transportation capacities are fixed. The case with uncertain demand, capacities, and costs is an important issue to research. In addition, the transportation distance can be utilized to identify the costs related to transportation. However, the type of transportation mode, vehicles traveling costs, and vehicle capacity can be used to define the transportation cost. Moreover, other metaheuristics such as NSGA-II , variable neighborhood search, particle swarm optimization, and ant colony optimization can be utilized to solve the proposed problem and to compare the developed algorithms. Furthermore, a non-dominated sorting technique could be incorporated into simulated annealing and tabu search to tackle the multi-objective issue. Another future work is to extend the supply chain models to address the parameters' uncertainty.

Practical implications and theoretical contribution
This study's theoretical contribution relates to the optimization of multi-objective supply chains. It integrated three objective functions: total supply chain profit, total risk, and total emission. The supply chain model includes several parameters and different types of decision variables. The presented model combined tactical/operational decisions and strategic decisions to provide a novel and more efficient supply chain system. This model also integrates several stages of the supply chain, which has a practical meaning in integrating supply chain performance and assessing supply chain partners. In addition, the proposed algorithms are designed to solve any large supply chain system and obtain near-optimal solution within acceptable computational time.
Practically, this work can be beneficial to design and configure many industrial supply chain networks. Some of these industries for which this work can be applicable include oil and gas supply chain, food supply chain, and petrochemical supply chain. The work in this paper helps managers to solve such complex supply chain systems and obtain a good configuration of the network. In addition, it helps the manager to act as a support system for the decision on the distribution of products across the supply chain network. Furthermore, it provides tradeoff plans from which decision-makers could choose based on their preferences.