A location-inventory model for distribution centers in a three-level supply chain under uncertainty

Article history: Received 26 August 2012 Received in revised format 14 September 2012 Accepted October 25 2012 Available online 28 October 2012 We study a location-inventory problem in a three level supply chain network under uncertainty, which leads to risk. The (r,Q) inventory control policy is applied for this problem. Besides, uncertainty exists in different parameters such as procurement, transportation costs, supply, demand and the capacity of different facilities (due to disaster, man-made events and etc). We present a robust optimization model, which concurrently specifies: locations of distribution centers to be opened, inventory control parameters (r,Q), and allocation of supply chain components. The model is formulated as a multi-objective mixed-integer nonlinear programming in order to minimize the expected total cost of such a supply chain network comprising location, procurement, transportation, holding, ordering, and shortage costs. Moreover, we develop an effective solution approach on the basis of multi-objective particle swarm optimization for solving the proposed model. Eventually, computational results of different examples of the problem and sensitivity analysis are exhibited to show the model and algorithm's feasibility and efficiency. © 2013 Growing Science Ltd. All rights reserved


Introduction
Supply chain consists of a network of suppliers, manufacturers, warehouses, distributors and customers who plan at changing raw material to final products, distributing them to customers and fulfilling their demand while minimizing (maximizing) cost (profit) of total chain.Commodities pass through different stages to be transferred from suppliers to customers, which may include different facilities (Min & Zhou 2002).In this regard, supply chain management plays an essential role in the cost reduction of companies and improvement of their competitive conditions.Supply chain network planning includes strategic, tactic and operational decisions.Meanwhile, strategic design of supply chain network, which is one of the most important elements of supply chain and can affect all of its decisions is very important and considerably affects chain planning and finally performance of companies.
One of the main decisions for designing a supply chain network is the issues related to location, which has a major part in investments.In a supply chain, a company seeks to locate facilities (such as plants, distribution centers and retailers) such that it can maximize (minimize) its total profit (cost).In the literature of this field, decisions of different strategic, tactical and operational levels are separately made due to the absence of accurate information about parameters of inventory cost and distribution at the time of making location decisions.For example, some researchers such as Zipkin (2000) studied and evaluated inventory policies in components of the chain by assuming that location decisions are clear or some researchers like Melo et al. (2009) developed retailers' locating models by ignoring the current tactical decisions.
But, the necessity of making high investments in these issues helped them convert to long-term (strategic) projects and, as mentioned above, they have long-term effects on the future costs and profits; as a result, it is very important that decisions of operational and tactical levels be considered while determining location decisions (Shen et al., 2003).In recent years, some researchers such as Shen and Qi (2007) showed that failure to consider the location and inventory costs at the same time of making decisions about location of facilities generally results in the generation of optimal sub-answers.As a result, the importance of making simultaneous decisions at different levels has led to development of integrated inventory -location models in recent years (Sourirajan et al., 2007) and many studies have been conducted to integrate tactical and strategic decisions and develop inventory-location models (Daskin et al., 2002;Shen et al., 2003).However, integration of these decisions at three levels has been less considered than the supply chain management attitude, majorly due to the consequences in operational levels and problems (Melo et al., 2009).
In addition, increase in competition in today's business environment, there are uncertainty associated with trade globalization and change in various factors of supply chain such as demand, supply and price.However, most of the available models assume fixed parameters and do not consider uncertainty.Change of location decisions is more difficult than that of inventory decisions which are more flexible.Therefore, location models should be able to include existing uncertainty in decision-making environment (Snyder, 2006).As mentioned, inventory-location models have been considered recently.One of the main assumptions of inventory-location decisions integration is the benefits resulting from risk pooling.In this case, inventories are kept in one place instead of storage in some different places, which causes reduction of inventory costs (Eppen 1979;Berman et al., 2012).Melo et al. (2009) performed a comprehensive review on the literature of supply chain management.Nozick and Turnquist (1998) studied how to include inventory costs (contingency reserve) in the classic problem of locating facilities of the plant and showed that these costs change almost linearly proportion to the number of distribution centers and, as a result, they can be regarded as fixed costs.Erlebacher and Meller (2000) had a broader view and studied a nonlinear integer inventory-location model by considering fixed costs, transportation and inventory costs.They presented a heuristic method to solve this class of problems using continuous space estimation and they applied it for 16 customer points.Shen et al. (2003) studied facility location problem in which facilities manage their inventory through policy of (r,Q).They developed location problem with fixed cost using EOQ approximation such that it included inventory costs.
The advantage of this model was that it benefits from integrative risk, which was mentioned by Eppen (1979) and the contingency reserve can be reduced to prevent from any potential shortage.This model is known as Location Model with Risk Pooling (LMRP) in literature and they used Column Generation method.Daskin, et al. (2002) presented an efficient solution based on Lagrangian relaxation approach, which solves the model at shorter time than the method of Shen et al. (2003) does.Teo and Shu (2004) studied the problem of logistic network considering inventory costs for multilevel locations of inventory storage.Of course, uncertainty of supply or demand was not been considered in their models.In another model, Romeijn et al. (2007) studied the previous model by adding term of contingency reserve (considering potential demand).Snyder et al. (2007) studied potential state of this model.In the literature of supply chain network design, there are few models, which have considered uncertainty in other parameters, except demand.Snyder et al. (2007) used scenario approach to include uncertainty of parameters.In fact, they formulated a potential planning problem to minimize the expected costs.
In addition, Qi et al. (2010) studied the effect of disorder in facilities for two-leveled supply chain, a supplier and several retailers, and sought to locate retailers optimally and allocate customers to them.In this model, disorder compensation periods followed exponential distribution and the lost sale was included.Chen et al. (2011) studied a reliable inventory -location model to optimize facility location decisions, allocation of customers and management of inventory in case distribution centers are at risk of disorder.They presented an integer planning model, which minimized total costs of their construction, expected value of customers and inventory holding costs under two normal scenarios and minimized confrontation with failure and used Lagrangian relaxation approach to solve the resulted model.Üster et al. (2008) studied single facility location, which includes location decisions, inventory completion (reordering periods), ordering costs, transportation and holding of inventory and presented three recursive heuristic methods as solution procedure.Park et al. (2010) presented a design model for three-level network by considering the contingency reserve.They also considered lead times depending on each pair of distribution center and supplier and location decisions of supplier/distribution center and inventory in an integrated way.Yazdani and shahanaghi (2010) presented a multi-objective probabilistic programming approach for locating distribution centers and allocating customers demands in supply chains, which considers risk in locating DCs, shipping products and also in arcs linking plant to DCs and DCs to customers though fuzzy parameters.The proposed model was solved by a probabilistic programming approach.Tancrez et al. (2012) considered inventory-location problem for three-level supply chain in which distribution decisions of distribution centers, allocation and amount of the transported commodity were made altogether and they were modeled as nonlinear continuous model.The proposed problem included costs of transportation, processing and holding of inventory.Berman et al. (2012) studied an inventory-location model in which distribution centers had inventory control policy (R,S) and found that its goal was to locate active distribution centers, allocate retailers to them and determine parameters relating to inventory policies while minimizing total costs.Tsao et al. (2012) studied the allocation-inventory-location integrated problem for designing a distribution network with several local distribution centers and retailers.Their modeling formulation located local distribution centers optimally and allocated retailers to them.The model also determined suitable inventory policy for each location while minimized total costs of the network and used Continuous Approximation (CA) and nonlinear planning for solving optimization problem.
In real supply chain systems, if the support facilities are considered for customers, reliability of supply chain and its general performance would be considerably improved.Rezaei et al. (2012) studied emergency response network design for hazardous materials transportation with uncertain demand, which is considered as fuzzy random parameter.They formulated the problem as a non-linear nonconvex mixed integer programming model.NSGA(ІІ) algorithm is applied to solve this model.Roghanian and Kamandanipour (2012) presented a closed-loop logistics network design based on reverse logistics models.A mixed integer linear programming is used to formulate this model.The problem considered single product and multi-stage logistics network for the new and return products and demand and rate of return are stated as fuzzy random parameters.The rest of this paper is arranged as fallows.In section 2, a brief description of robust optimization is presented.In section 3, the given problem is formulated as a robust optimization model.In section 4, a solution approach based on meta-heuristic algorithm is explained.Section 5 presents the numerical examples and sensitivity analysis for the most important parameters.Conclusion and future research of this study is provided in section 6.

Robust optimization
Despite the fact that there are many applications of stochastic programming, routine stochastic programming models acutely are restricted because of disability in managing risk aversion or decision makers' preferences directly.Here, we used a robust programming, introduced by Mulvey et al. (1995), which is an improved stochastic programming to deal with the preferred risk aversion of decision makers, which was not possible to use in routine stochastic programming (Bozorgi Amiri et al., 2012;Azaron et al., 2008) to locate distribution facilities in a three levels supply chain network under uncertainty.
In this method, the variability term was supplemented to the main objective function by a related weighting parameter in order to show the tolerance of modelers' risk.In the remaining, a concise description of robust optimization is presented.Let x be first stage (design) variable vector and y s be the second stage control variable vector.Let M, N, O be matrices of parameter and p, q be vectors of parameters.Let M, p be certain parameters and N, O, q be uncertain ones.Let S be the set of scenarios to model uncertainty with associated probability of occurrence of each scenario (p s ): S={1,2,…,s} and ∑p s =1 so we have N s , O s , q s .The possible infeasibility of model is presented by δ s .(ifmodel is feasible 0, otherwise a positive value is calculated by Eq. ( 3) model formulation is as follows (Bozorgi Amiri et al. 2012): 1 s 1 min (x,y ,...,y )+ p( ,..., ) The first term in Eq. ( 1) shows the solution robustness, which seeks for less cost and risk aversion degree.The second term indicates the model robustness, which gives penalty solutions by unmet demands or exceeding from each physical constraint.ω is a weight measuring the trades-off between the first and the second term of Eq. ( 1).According to Bozorgi amiri et al.'s work, we use υ s =z(x,y s ) as the main objective function under scenario s.The solution is a high-risk decision when the variance of υ s is high.Mulvey et al. (1995) used a quadratic form of variance, which is nonlinear and complicated.
To handle this difficulty, here we use an absolute deviation as Yu and Li (2000) proposed, which is as follows, ( ) where λ is the weight of the less sensitive-solution to data changing in all scenarios.For minimizing the Eq. ( 5), Yu and Li (2000) presented an effective method, which is modeled as follows, min 2 If υ s is bigger than ∑p s υ s , θ s is equal to 0, otherwise θ s =∑p s υ s -υ s .In this study, We use Yu and Li's method but as the expected value of costs and their variance are in contrast, we form a two-objective model, which separates the two presented terms in Eq. ( 6) in order to enhance the model efficiency and find Pareto solutions.This provides a good condition for decision makers to make decision according to their preference.

Problem description
In this section, an integrated two-objective robust inventory-location model is presented to design three-level supply chain network including suppliers, distribution centers and customers as integer nonlinear programming in case of uncertainty, which integrates location, inventory and allocation decisions.The presented model's assumptions are as follows:

Assumptions
1-More than one product can be supplied in the chain, each one of the products has different volume, procurement, shortage, holding, ordering, and transportation cost.
2-Capacity of the suppliers and distribution centers is subject to uncertainty due to events such as fire, earthquake, etc. 3-There is uncertainty in parameters such as demand, supply, purchase price and costs and some discrete scenarios belonging to the set of possible scenarios S are used for showing this problem.4-There are some candidate points for establishment of distribution centers, which are selected for each location's fixed cost.5-Each distribution center can be constructed only in one of the available sizes (small, medium and large).6-Each distribution center can be supplied from suppliers and other distribution centers (logistic cover) if possible.7-There is no limitation of single-facility supply for customers and distribution centers.8-Inventory is kept only in distribution centers.In this case, this inventory is fined.9-Remaining inventory holding cost is different depending on whether the remaining inventory is related to supplier or logistic distribution center.10-The package shortage cost differs depending on whether shortage results from failure to select customer for providing services because there is uneconomical or insufficient distribution center to provide services while the model is as the lost sale.11-In this model, there is inventory ordering policy (r,Q) and EOQ approximation approach has been used to determine its parameters based on the work of Axsater (2006).It has also risk pooling property.
The goal of this model is to design a distribution network, which is solved to specify location and number of distribution centers, inventory order amount in each one of the distribution centers, the allocation of customers to distribution centers as well as distribution centers to suppliers by aiming at minimizing the expected value of costs and variance of these costs.Variables of the first stage (design) and the related fixed location cost were final but variables of the second stage (control) and its related parameters such as demand, supply, etc. are assumed uncertain.Uncertainty is captured by some specific discrete scenarios.In the remaining, the symbols relating to this problem are presented.h̒ jc : holding cost of a commodity unit c in distribution center j (in case it receives commodity from supporting distribution centers in addition to suppliers) α : in case the support is not used, it equals 0 and in case the support is used, it equals 1 M: a large number

Nondeterministic Parameters
q js : a percent of capacity j which remains active under scenario s q is : a percent of capacity i which remains active under scenario s φ ics : purchase cost of a commodity unit c from supplier i under scenario s φ j'cs : purchase cost of a commodity unit c from supplier j under scenario s n jcs : number of orders for commodity c by distribution center j under scenario s s θ : the variable applied for linearization of absolute deviation of costs

Mathematical Formulation
In this section, a new robust mathematical model is presented in which uncertainty is expressed using a finite number of discrete scenarios.As mentioned before, EOQ approximation approach was applied here to use policy of (r, Q) and determine the number of order of distribution center per year.
On the other hand, considering that shortage was permissible but irrecoverable, it was proved that, in the presence of inventory system in this state, economic order was calculated through Wilson relation.
Reordering point was obtained as follows, considering that demand was specified in each scenario: (10) LT= LT Of course, since the number of ordering was optimized in the model presented in this thesis, economic order can be calculated through the following formula by having optimal number of ordering: In this regard, ordering policy parameters (r,Q) can be calculated.Then, mathematical model of the problem and a scheme of inflows and outflows are presented for one node of distribution centre j, according to Fig. 2 in order to understand the problem.
Eq. ( 12) is the first objective function, which includes minimization of fixed location costs and expected value of costs of purchase, transportation, inventory holding, shortage in distribution centers, lost sale and ordering along the one-year planning horizon.Eq. ( 13) is the second objective function of this problem, which is variance of purchase, transportation, inventory holding, shortage and ordering costs.This has been considered as absolute magnitude and it is linearized based on the available literature in this field, like what was referred to in the study by Yu and Li ( 2000).Constraints of the problem were as follows: Subject to: ) 14 ( j,c,s ijcs j jcs jl jkcs jj cs j l jcs jkcs Eq. ( 14) is inventory balance equation in distribution center for any kind of commodity.Eq. ( 15) is included to calculate the number of ordering along the planning horizon.Constraint ( 16) shows that distribution center j provides services to other distribution centers as supporting center when it has been established with size of l in location j.Constraint ( 17) also reveals that other distribution centers (j̒ ) provide service to distribution center j as a logistics center when location j of distribution center with size l is constructed.Constraints (18) and constraint (19) show dependency of supplier i and customer k on distribution j.Constraint (20) depicts that total volume of the commodities, which distribution center j can deliver to customers and other distribution centers is the same as its accessible (active) capacity.Constraint ( 21) is constraint of capacity of supplier i considering its active capacity.Constraint (22) indicates that, at most, one with any possible size is constructed for the distribution center.Constraint (23) shows that if the amount of the commodity suppliers deliver to distribution center j is lower than its demand, it will receive service from logistic distribution centers shown using binary parameter of α.Constraint ( 24) is a balance equation for node k (customer).Constraint ( 25) is constraint of the problem variables and constraints ( 26) and ( 27) are the constraints resulting from linearization of costs variance.

Solution procedure
In this paper, first, the model is solved using Lingo 9 software and using Epsilon Constraint Method; then, MOPSO meta-heuristic algorithm is used to solve the resulted problem for larger problems.Each one of them is briefly explained below.

Epsilon Constraint Method
Epsilon Constraint Method is one of the well-known approaches for handling multi-objective problems, which solve such problems by transferring all objective functions into constraints and keeps only one of them in each phase as objective function (Ehrgott & Gandibleux, 2002).In this case, Pareto Border can be created with ε constraint method (Bérubé et al. 2009).
The following summarizes the necessary steps of ε -constraint method, 1-One of the objective functions is selected as the main objective function.
2-Considering one of the objective functions, the problem is solved and optimal values of each objective function are obtained considering one of the objective functions.3-The interval between two optimal values of objective sub-functions is divided into the predetermined number and a schedule of values is obtained for ε 2 ,..., ε n. 4-The problem is solved at any time with main objective function with each value of ε 2 ,..., ε n. 5-The found Pareto's answer is reported.
In order to study the reason for considering this model as bi-objective, problem 3 is considered with the problem codes of 2-1-3-2-3.Based on steps of the problem, 10 sub-problems are generated in which optimal value of z 2 is determined at any time by selecting the second objective function as the main objective function and putting the first objective function in constraint based on ε 2 ,..., ε n..The results obtained from these calculations are given in Table 2.

MOPSO Algorithm
Considering the fact that the proposed model can be solved in the simplest state as an allocationlocation form without capacity constraint, which is NP-hard based on the work of Megiddo and Supowit (1984), this model, which is the development of the mentioned basic model is NP-hard.For this reason, meta-heuristic algorithm of multi-objective particle swarm optimization (MOPSO) with had high convergence speed is used.

Introduction of PSO Algorithm
PSO is one of the population based optimization algorithms, which was presented by Kennedy and Eberhart (1995).First, the set of particles is placed in the response space and starts moving with initial velocity.Then, these particles move in response space and are evaluated according to special criteria in each stage.Over time, these particles accelerate with specified velocity towards other particles available in their communication group in multidimensional search space, which have higher fitness value.Any position of particles shows a solution for the problem.Notion of the algorithm are introduced in Table 2.

Table 2
Notations of PSO algorithm G :fitness value function V max :maximum velocity that each particle can take.if V> V max THEN V= V max .
v id (t) : velocity of particle i in d th dimension in t th stage P id :location of the best position found so far by each particle in dimension d ( ) x t : location of each particle vector P gd :: location of the best position found by all particles in dimension d W (INTERIA WIEGHT) : controls intensification and diversification L P : vector of best found position of each particle so far d : dimensions of solution space g P r : vector of best found position of all particles x id (t) : location of particle i in d th dimension in t th stage C 1 and C 2 : fixed coefficient to control the impact of P id and P gd ( ) P t r : vector of best found solution 1 2 and φ φ : uniform random number between [0,1] At any iteration, velocity and position of particles are updated according to Eq. ( 29) and Eq. ( 30): )

Introduction of MOPSO Algorithm
Traditional PSO cannot be used for multi-objective optimization problems and considering that a set of solution are presented instead of one in the multi-objective models, changes in algorithm are needed so that one can reach this set of solution called Pareto Optimum.Generally, three objectives should be considered in order to solve multi-objective problems: (Eckart et al. 2000).
1-Increasing the number of found Pareto points, which is called quality metric.2-Decreasing the resulted Pareto border distance by the algorithm with global optimal Pareto border called spacing metric.3-Increasing the range of the obtained solutions such that there is a uniform distribution of the response vectors as much as possible called diversity metric.
Fig. 2 shows proposed multi-objective PSO algorithm pseudo code of this paper.In fact, binary MOPSO method of zero and one was used.Since space of the problem was binary, binary algorithm which was presented by Eberhart and Kennedy was used in this research (Kennedy & Eberhart 2001).
In order to update archive of Pareto solution, Roulette wheel operator was used.Most of evolutionary meta-heuristic algorithms use a random approach for producing initial solution.Here, this approach was followed to produce initial solution (Alvarez-Benitez et al. 2005).

parameters tuning of the Applied MOPSO Algorithm
In the meta-heuristic algorithm of multi-objective particle swarm optimization, parameters of the algorithm are set using trial and error method after consecutive periods of algorithm execution as well as considering the results obtained from the related papers in this field.For determining the number of initial population, ϕ 1 and ϕ 2 , and the maximum iteration which is stopping criteria, response surface method (RSM) in design expert software was used.
In MOPSO, initial population is randomly produced with uniform distribution.Primary population is considered 60 in each repetition.The maximum number of repetitions at any time of the algorithm execution is considered 60.Values of ϕ 1 and ϕ 2 are both considered as 2.05.Personal and ultimate learning coefficients of c 1 and c 2 are as a fixed coefficient (k=l) of ϕ 1 and ϕ 2, respectively.w damp is considered 0.999 and α is considered as 0.1.Here, stop condition is achieved when algorithm reaches the maximum number of repetitions.

Numerical Examples
In this section, first, lingo 9 was used in order to validate the model.Then, since this software is not able to solve multi-objective problems, the conflict between two objectives of the problem related model is revealed by using Epsilon Constraint method, which is described in Section 4.1 to specify the reason for considering the bi-objective model.Since the model is NP-hard, meta-heuristic algorithm of multi-objective particle swarm optimization (MOPSO), which had high convergence speed has been used.Then, 10 designed problems were solved using both exact and meta-heuristic methods and its calculation results have been comparatively shown.Lingo9 software was used for the exact method and MATLAB 2009 software was used for meta-heuristic method; all the tests were performed by Pentium3 (Intel (R) Core TM 2CPUT5500@1.66GHz)and Microsoft Windows XP Professional operating system.Then, specifications relating to these 10 problems are given in Table 3.

Studying Conflict between objective Functions
Here, problem 3 is selected in order to study the conflict between objective functions and has been solved based on what was referred to Section 4.1 using Epsilon Constraint method; the results are given in Table 4. Table 4 shows that the assumed objective functions, of which one includes minimization of the expected value of supply costs and another one includes minimization of variance of the supply costs, change in the opposite directions; when one is improved, another one will be worsen.As a result, these two objective functions are conflicting and putting them beside each other does not seem reasonable.While separation of these two objective functions causes the manager to have more focus on the variance of costs, which is mostly neglected in decision making.

5.1.The Results Obtained from Numerical Examples
In this section, 10 problems are randomly designed (considering information available in Table 3) and solved with both exact and meta-heuristic methods mentioned in Section 4 to study efficiency of these methods.The obtained results are given comparatively in In case the objective functions are improved, the minimum and maximum percent of improvement are calculated from 1% to 11% for the first objective function (Z 1 ).In addition, the minimum and the maximum percentage of improvement are from 12% to 21% for the second objective function (Z 2 ), respectively, which is regarded as considerable improvement compared with problem 7 with the worsened objective function.In the problems with large size, it is found that exact solution in lingo is not able to find the response while acceptable solution are created using MOPSO meta-heuristic algorithm.
Both exact method and MOPSO meta-heuristic algorithm have reached reasonable solutions in a short period of time for small instances.For large scale problems, no response is obtained within 48 hours time execution using Lingo software.However, MOPSO meta-heuristic algorithm has been able to create justifiable solution within reasonable amount of time.Considering the mentioned cases, one can conclude that this algorithm is more efficient than the presented exact solution and is able to present justifiable and acceptable answers for large sizes.

Sensitivity Analysis
In order to study the effects of some important parameters on values of objective functions which are determined by trial and error method, five parameters including fixed cost of location, capacity of distribution centers, capacity of suppliers, shortage cost and demand rate have been selected and executed in MATLAB for problem 3 and their results are given in Fig. 4. As it is depicted, spacing metric has a non-decreasing trend with the increase of values of the intended parameters for shortage cost and the capacity of distribution center has an ascending trend for the fixed cost of location and has a fluctuating trend in other cases.In diversity metric, it has an ascending trend for the fixed cost of location and demand rate and has a descending trend for capacity of supplier and distribution centre and shortage cost and this trend holds true for quality metric.Regarding mean of the first objective function, an ascending trend is found in the fixed cost of location, shortage cost and demand rate and a descending trend is found in capacity of supplier and distribution center.As far as the second objective function is concerned, there is a generally descending trend in the fixed cost of location, capacity of supplier and distribution center and there is an ascending trend in shortage cost and demand rate.At the end, one can conclude that the related parameters influence the values of objective function and comparative indices of MOPSO, which are determined by considering the conditions related to the given supply chain and utility of the decision makers.However, it is worth noting that values of objective function changed a little in all the cases despite the changes of parameters and, as shown in the figures, slope of the diagrams relating to mean of the objective functions is mild, which indicates robustness of the applied algorithm.

Conclusion and future research
This paper has concentrated on the joint location-inventory problem where a three-level supply chain networks under uncertain environment through finite discrete scenarios.We proposed a robust optimization model and formulated the model as a two-objective mixed-integer nonlinear programming.In addition, a multi-objective particle swarm has been applied to solve the model.Then computational trials for three sizes of the problem (small, medium, large) and MOPSO algorithm have been compared with the mean of obtained solution of lingo for each objective function by the εconstraint method.Numerical instances pointed out the MOPSO algorithm performed more logical compared with other algorithms.Meanwhile, MOPSO could find a near-optimal solution for large instances, where optimal solutions are not available.All in all cases, it could be said that the proposed algorithm has been able to solve the problem in reasonable amount of time, efficiently.Finally, a sensitivity analysis has been presented to survey the impacts of the most important factors on the objective functions and specific metrics of MOPSO algorithm, which indicates based on the MOPSO, the two mentioned objectives, show slight changes that approve the robustness of applied algorithm.
The main contributions of this research lie in integrating strategic (location) decisions with tactical ones (inventory and support system).We have presented a robust optimization model by considering different uncertainties such as uncertain demand, supply, procurement, transportation cost, orders' amounts and facilities' capacities.This helps us design more realistic and flexible supply chain networks, contemplating the variance of total cost as an objective function.
At last, as the future research directions of this work can be as follows: (1) merging the proposed model by pricing policies.(2) developing the model for more than three level supply chains and considering other uncertain factors such as lead times in order to more compatibility by real-world's situations.(3) taking into account the discounting policies about procurement, transportation cost, and etc. (4) using other function such as financial risk or reliability as an individual objective.

I
: set of suppliers J: set of distribution centers K: set of customers L: set of the assumed sizes for distribution center (small, large, medium) S: set of possible scenarios C: set of all demanded commodity 3.3.Deterministic Parameters f jl : fixed cost of opening distribution center j with size l F jc : fixed ordering cost of any distribution center j for each commodity c P s : probability of occurrence of scenario s π kc : shortage cost in distribution centers for one unit of commodity C resulting from demand of customer k (the penalty of un-met demand of assigned customer to a DC due to uncertainty of supplier) π̒ kc : shortage cost resulting from failure to allocate some of the customer's demand k for each commodity unit C (lost sale) υ c : per unit required space for each commodity C S ic : amount of commodity c which supplier i can supply Cap l : different type of opened DCs' capacities (in cubic meters) h jc : holding cost of a commodity unit c in distribution center j (in case it receives commodity only from suppliers) C ijcs : transportation cost from supplier i to distributor j for each commodity unit c under scenario s C j'jcs : transportation cost from logistic distributor j' to distributor j for each commodity unit c under scenario s C jkcs : transportation cost from distributor j to customer k for each commodity unit c under scenario s d kcs : demand of customer k for commodity c under scenario s3.5.Continuous and Binary VariablesX ijcs : amount of commodity c transported from supplier i to distribution center j under scenario s Y jkcs : amount of commodity c transported from distribution center j to customer k under scenario s I jcs : inventory of commodity c which is stored in distribution center j under scenario s b jkcs : shortage of commodity c resulting from demand of customer k in distribution center j under scenario s b̒ kcs : shortage of commodity c for customer k under scenario s which has not been allocated to any distribution center because of limitation and uncertainty of capacity of distribution center Z jl : it is 1 in case distribution center is opened with size L in location j; otherwise, it is 0.

Fig. 2 .
Fig. 2. Input and output flows of DC

Table 1
Characteristics of reviewed articles about facility location under uncertainty

Table 3
Data of generating parameters of numerical examples

Table 4
Results of applying ε-constraint for test problem no. 3 in lingo 9

Table 5 .
As we can observe from the results of Table5, for small instances, global solutions are available but for large scale problems, we need to use MOPSO.Using Eq. (34), we may compare the gap between the performance of B&B problem with MOPSO.The gap varies from 2% to 12% depending the size of the problem.