A novel robust chance constrained possibilistic programming model for disaster relief logistics under uncertainty

Article history: Received November 4 2015 Received in Revised Format December 21 2015 Accepted February 25 2016 Available online February 25 2016 In this paper, a novel multi-objective robust possibilistic programming model is proposed, which simultaneously considers maximizing the distributive justice in relief distribution, minimizing the risk of relief distribution, and minimizing the total logistics costs. To effectively cope with the uncertainties of the after-disaster environment, the uncertain parameters of the proposed model are considered in the form of fuzzy trapezoidal numbers. The proposed model not only considers relief commodities priority and demand points priority in relief distribution, but also considers the difference between the pre-disaster and post-disaster supply abilities of the suppliers. In order to solve the proposed model, the LP-metric and the improved augmented εconstraint methods are used. Second, a set of test problems are designed to evaluate the effectiveness of the proposed robust model against its equivalent deterministic form, which reveales the capabilities of the robust model. Finally, to illustrate the performance of the proposed robust model, a seismic region of northwestern Iran (East Azerbaijan) is selected as a case study to model its relief logistics in the face of future earthquakes. This investigation indicates the usefulness of the proposed model in the field of crisis. © 2016 Growing Science Ltd. All rights reserved


Introduction
A disaster is an unscheduled, overwhelming incident in association of people with their environment causing death, injury, and extensive casualties (Rubin et al., 1985).Disasters can be categorized as natural and man-made.Examples of the first one includes earthquake, flood, storm, drought, hurricanes and terrorist attacks, chemical leakages are some instances of man-made (Caunhye et al. 2012).From 2000 to 2007, the number of reported natural disasters around the world was approximately 460 disasters per year, which cost between 100 million and 400 million victims per year (Haghani & Afshar, 2009).Although natural disasters are unexpected, their damages can be minimized with proper preventive plans.Relief distribution center (RDC) location and relief distribution are among important strategies to improve the relief performance, since the numbers and the locations of RDCs, and the amount of supplies pre-positioned in them will directly affect the response time and the cost of logistics.So, this creates motivation to model the RDC location, and the inventory decisions associated with the relief distribution in a relief network.One of the important features of a relief chain, which could increase the complexity of the disaster relief logistics, is associated with the uncertain and dynamic factors existing in the postdisaster environment placed into three major groups (Bozorgi-Amiri et al., 2013): 1.The uncertainty in relief supply caused by the delays in the supplier deliveries, unknown usable resources due to the roads and the building's destruction, and the unpredictable involvement and contribution of the suppliers, 2.The uncertainty in relief demand; mainly due to the inaccurate assessments or the volatility of the demand, due to the people self-sufficiency improvement, people movements to gain more relief aids, and disease outbreak.3.Uncertainty in relief costs; resulted as the uncertainty associated with routs, suppliers, etc.
To cope with all the above mentioned uncertainties, this paper applies a possibilistic programming method to model the relief distribution system; in which all the uncertain parameters are considered in the form of fuzzy trapezoidal numbers.This model tackles the problem as a multi-objective, possibilistic, mixed-integer, nonlinear programming model.Then the problem is then solved using the LP-metric and the Improved augmented -constraint (AUGMECON2) methods.
The main contributions of this paper can be presented as follows: • The proposed method uses a robust chance constrained possibilistic programming to cope with the uncertainties of a disaster relief logistics network.
• It provides a three-objective mathematical model which simultaneously takes into account the distributive justice in relief distribution, the risk of relief distribution, and the total logistics costs.
• It considers the relief commodities priority and affected areas (AAs) priority simultaneously.
• It takes into account not only the probability of facility (suppliers and RDCs) disruption during the disaster, but also proposes a new work by considering the disaster retrofitting of distribution centers buildings (in preparedness phase).
• It looks for the uncertainty of relief supply in a new fashion, by distinguishing between the suppliers' supply ability in preparedness and response phases.So that the pre-disaster supply amount can be abundant and deterministic, while the post-disaster supply amount is usually limited and uncertain (due to the emergency situations).
• Finally, it exerts the distribution standard of the relief organization throughout the network, by meeting a minimum percentage of each affected area's demand level.
The rest of this paper is structured as follows: the relevant literature is reviewed in Section 2. In Section 3, the concept of the robust programming is proposed.Then, problem statement, notation and mathematical model are given in Section 4 followed by the description of the solution methods provided in Section 5. Evaluation of the proposed robust model is provided in Section 6.In addition, introduction of the case study and its experimental results are provided in Section 7. Finally, concluding remarks are stated in Section 8.

Literature review
Facility location literature is very broad and rich topic since it considers strategic decisions for a wide range of public and private plants.The models in this area can be placed into four main groups (Owen& Daskin,1998): Deterministic location problems, Dynamic location problems, Stochastic location problems, which consists of Probabilistic models and Scenario based models, Fuzzy location problems, that are grouped into Flexible programming and Possibilistic programming problems.
Decision making based on a deterministic model may increase the existing risks and can make the tough emergency situations of a disaster relief even more disastrous.Among the above four categories, only the stochastic and the fuzzy programming methods can include the uncertainty of the parameters into the model.
There are several review articles in facility location literature, written in the past few years (Caunhye et al., 2012;Luis et al.,2012;Kovács& Spens,2007); providing a comprehensive overview of the humanitarian relief logistics models.Most models in the facility location literature combine the problem of relief facility location (establishing a new facility or choosing from existing facilities), with relief commodities pre-positioning, evacuation, and relief distribution problems (Jia et al., 2007).Also, most models in the literature of relief facility location are based on mixed integer programming with binary location variables.Furthermore, as the relief facility location models are all used for predisaster planning, they are all found to be single-period (Caunhye et al., 2012).In some of the recent studies, researchers have motivated to address stochastic optimization in relief facility location planning involving facility locating and distribution of emergency commodities by probabilistic scenarios representing disasters and their outcomes (see Table 1).On the other hand, four major disadvantages are identified for the stochastic programming approach (Ben-Tal & Nemirovski, 2008) : A. In most cases due to the lack of sufficient historical data, the actual and precise distribution function of uncertain parameters cannot be found.B.The large number of scenarios used to exemplify the uncertainty of the relief environment, may contribute to the computational complexity of the problem (Caunhye et al., 2012).C.This approach cannot take decision makers risk-averse behavior directly into the model, so despite using this method in a large group of existing models they have a limited application.D. On stochastic optimization, scenarios are formed based on possible deterministic observations of the uncertain parameters and the solution is generated based on those scenarios.So the answer might become infeasible due to other observations of the uncertain parameters.Observations that although their occurrence probability is really low, but their occurrence will impose a high price to the entire relief network.
Also, a risk-averse decision-making approach that is able to amend the third objection of the stochastic models, is the Robust optimization theory.This approach, which is first presented by Meloy et al. (1995), can be applied in both stochastic and fuzzy programming methods in the face of uncertainties to exert the decision maker's risk aversion attitudes into the modeling.In recent decades, many researches have addressed location problems by applying fuzzy logic methods.For example, Bhattacharya et al. (1992) used a fuzzy goal programming method to solve their model.Canos et al. (1992), Darzentas (1987), Rao & Saraswati (1988), all addressed fuzzy location problems, but they all considered deterministic parameters for their models.Also Zhou & Liu (2007) located facilities and allocated demand points to them considering fuzzy demands and facility capacity constraints.
In this paper, we consider the disadvantages of stochastic programming and present the application of fuzzy theory in representing the uncertainties of relief environment.To the best of our knowledge, it is the first time in relief facility location literature, that the chance constrained possibilistic programming method is applied and uncertainty of supply, demand and the costs of the relief environment are considered in the form of fuzzy trapezoidal coefficients.Then the robust optimization approach is applied to involve the decision maker's risk aversion attitude in our model.

Robust programming (RP)
Robust programming provides risk aversion methods in dealing with uncertainty in optimization problems.A solution of an optimization problem is called a robust solution if it simultaneously fulfills two types of robustness: "solution robustness" (the solution is nearly optimal for all possible values of the uncertain parameters) and "model robustness" (the solution is nearly feasible for all possible values of the uncertain parameters) (Mulvey et al.,1995;Ben-Tal& Nemirovski,2002).

Robust possibilistic programming (RPP)
Robust possibilistic programming is a novel possibilistic approach provided by Pishvaee et al. (2012), which utilizes the possibilistic programming and fuzzy logic concepts of (Inuiguchi & Ramık, 2000;Liu & Iwamura, 1998;Dubois & Prade, 1997;Heilpern, 1992), and integrates them with the robust programming frameworks.To illustrate the RPP approach, first consider the following typical singleobjective fuzzy model: where the vectors f, c, d, and the matrix N represent uncertain parameters in terms of fuzzy trapezoidal numbers and vectors y and x denote the binary and the continuous variables of the model, respectively.
To make the basic form of the chance constrained programming (CCP), Pishvaee et al. (2012) used the expected value operator in the objective function and the necessity measure in the chance constraints which include imprecise parameters.Thus, the basic form of the CCP model can be stated as follows: According to Pishvaee et al. (2012) the equivalent deterministic form of the above formulation is as follows: where the first term in the objective function is the expected value of the objective function of model 1, which is the minimization of the expected total performance of the concerned system.The second term, i.e. , denotes the difference between the two extreme possible values of z.This term controls the solution robustness of the solution vector by minimizing the maximum deviation over and under the expected optimal value of z.In other words, is associated with the weight (importance) of this term against the two other terms of the objective function and and are defined as follows: The third term in the objective function of model 3, i.e., 1 , determines the confidence level of the first chance constraint (first constraint on the model 3).Also, is the penalty unit of possible violation of the first chance constraint, and 1 , is the difference between the worst case value of this constraint (based on the range of its imprecise parameter) and the value that is used in this chance constraint.In fact, the third term controls the model robustness of the solution vector.Similarly, the forth term, i.e., 1 , determines the confidence level of the second chance constraint (third constraint of the model 3) and controls the model robustness of the solution vector.Also, in model 3, variables , ∈ 0,1 are the confidence levels of the chance constraints and Pishvaee et al. assumed , 0.5 to satisfy the chance constraints with a chance level greater than 0.5.
Finally, as the model 3 is a non-linear model, it is converted to an equivalent linear form as follows (see Pishvaee et al. (2012) for more details): There are some cases where the decision maker is not sensitive to both over and under deviations from the expected optimal value of the objective function, therefore, Pishvaee et al. ( 2012) changed the second term in the objective function of the above-mentioned model to form a new model as follows: where, F is the feasible region of the model 6.

Problem statement
We consider a three-stage disaster relief logistics network in which the first stage contains suppliers, the second one consists of relief distribution centers (RDCs) and the last one is associated with the set of affected areas (AAs).As a pre-disaster planning, we assume locating RDCs from a set of known candidate locations so that their storage capacity and disaster retrofitting decisions are also determined by this selection.Also, in the preparedness phase, frequently used relief for commodities are pre-positioned in the selected RDCs to accelerate the after-disaster relief operations.In the response phase, we consider the commodity transportation from suppliers to RDCs, between RDCs (backup coverage) and from RDCs to the AAs.We make the following assumptions to model this problem: (1) The location of candidate RDC points and potential AA points are identified by the decision makers before the planning time.
(2) The capability of suppliers and RDCs might be partially disrupted by a disaster (3) The uncertainty of supply, demand and the costs of the relief environment are considered in terms of trapezoidal fuzzy coefficients.(4) Three types of relief commodities (water/ food/ shelter) are supposed to be delivered so that each type has its own volume and cost of procurement, storage, and transportation.( 5) An RDC can be opened with only one of the three possible storage capacities, small, medium, large, and seismic retrofitting levels (not retrofitted, partly retrofitted, totally retrofitted); subject to the associated setup cost.(6) To ease the relief coordinations, each AA only serves with one RDC.(7) In the response phase, not only the commodity shortages in AAs, but also the excess inventory stored in RDCs, are penalized.

4.1.Mathematical model
The According to the robust chance constrained approach described in the previous section, the proposed model has the following form: First objective: The first objective function is associated with increasing the distributive justice in relief distribution.In this case, it minimizes the summation of weighted maximum relative lack of any kind of relief commodities in AAs, by considering the priority of each commodity type and the priority of each AA in relief distribution.It is the first time in disaster relief logistics, that the distributive justice in relief distribution is considered along with commodities' priority and AA's priority in getting relief services.The linear form of this objective function in the form of chance constrained programming is as follows: As the Eq. ( 10) is a chance constraint, to provide model robustness of the solution vector, will be added to the final robust objective function: Second objective is also as follows, The second objective minimizes the risk of relief distribution and implicitly increasing the relief distribution speed from RDCs to AAs.The linear form of this objective function is reported as the form of Eq. ( 14) to Eq. ( 18): .
Third objective: The third objective of the problem minimizes before and after disaster logistics costs.Eq. ( 20) consists of construction and retrofitting costs of the RDCs, procurement costs of the relief commodities and transportation cost of pre-positioning the relief commodities in RDCs in preparedness phase: Also, Eq. ( 21) encompasses after disaster costs, such as procurement and transportation costs of the relief commodities throughout the network in response phase.According to the adopted robust chance constrained approach, the deterministic form of the third objective function is as follows: ., This constraint controls the commodity balance of each RDC.In this constraint, the term cjre pj is obtained from the following equation: In order to linearize the Eq. ( 25), we define variables 1 2 3 , , cj cj cj w w w as follows: So that we have: Accordingly the linear form of the first constraint is equivalent to Eqs. (30-41): , , (41) Also, Eq. ( 31) ensures that at most one RDC can be constructed at each RDC candidate point.Constraint 2: . , This constraint shows the relief organization standard in responding to the affected areas of demands and Eq. ( 43) shows its deterministic form: (43) Afterwards, 2 CLCC will be added to the final robust objective function to fulfill model robustness of the solution vector: Eq. ( 45) to Eq. ( 47) ensure having backup coverage between only any two established RDCs.Also, Eq. ( 48) prevents suppliers from transferring commodities to a not opened RDC and Eq. ( 49) guarantees sending commodities to AAs only from an opened RDC.
These two constraints ensure that commodity pre-positioning and storing additional commodities in RDCs must be based on their storing capacities.Constraint 10: 1 , According to this constraint, in preparedness phase, the amount of commodity c procured from supplier i, cannot exceed the supplier's capacity.Constraint 11: . 2 , This constraint ensures that, in response phase, the dispatched commodity from each supplier is limited by its usable inventory amount.Also, the deterministic form of this constraint, in chance constrained programming, is as follows: Moreover, in order to fulfill the model robustness of the solution vector, the term 3 CLCC should be added to the final robust objective function: The commodity inventory of the suppliers in preparedness phase is supposed to be a deterministic parameter, Since before disaster they have enough time to procure and supply an abundant amount of commodities.While in response phase, as suppliers must emergently supply relief commodities, their commodity inventory is supposed to be a non-deterministic parameter.So suppliers' capabilities in procuring commodities are different in preparedness and response phases, while to the best of our knowledge, this note is ignored in all the previous studies in the disaster relief logistics literature.So that, in all the relief logistics models considering the uncertainty of the relief supply, the relief supply in preparedness phase and response phase are both supposed to be an uncertain parameters.Constraints 12-14: (56) 6 . , According to the above three constraints, each AA allocates to exactly one RDC and receives reliefs only from its dedicated RDC.This allocation can also ease the relief coordinations while having backup coverage between the RDCs of the relief network.

Solution methods
First, Lp-metric method is utilized, under which the proposed model becomes a single-objective problem and it enables us to evaluate the proposed robust model against its equivalent deterministic form.Then, in the next step, to provide a comprehensive sight of our three-objective problem, augmented -constraint method (AUGMECON2) is presented and this method is used to solve the case study.

LP-metric
In this method, first, the optimal value of each objective function is calculated by solving the relevant one-objective problems.Then, a single objective function is employed to minimize the summation of normalized differences between each objective function and its optimal value, on the solution space of the multi-objective function (Soltani et al.,2015).
We used this method by minimizing the following objective function on the solution space of our model: ])

Improved augmented -constraint method (AUGMECON2)
The -constraint method solves a multi-objective problem by optimizing one of the objective functions while using the other objective functions as constraints, incorporating them in the constraint part of the model.In order to apply this method, it is necessary to have the range of the objective functions used as constraints.Payoff table is a common approach to calculate these ranges.This table is made with the result of individual optimization of the objective functions while the nadir value is usually approximated with the minimum of the corresponding column (Steuer, 1986;Miettinen, 2012).
For the implementation of the ordinary -constraint method two points must be considered: first, the range of the objective functions over the efficient set, mainly the calculation of nadir values.Second, the ,  guarantee of efficiency of the obtained solution (Mavrotas, 2009).To overcome these ambiguities, augmented -constraint method (AUGMECON) is proposed (Mavrotas, 2013).In AUGMECON method, lexicographic optimization is used to construct the payoff table with only Pareto-optimal solutions.There is a trade off between the density of the produced efficient set and the computation time, so we can control the density of the efficient set.In this method, in order to guarantee the efficiency of the obtained solution, the objective-function constraints are transformed to equalities by incorporating the appropriate slack or surplus variables.These slack or surplus variables are used as the second term, with lower priority in a lexicographic manner, in the objective function, forcing the program to produce only efficient solutions (Mavrotas, 2013).According to (Mavrotas, 2009), in the AUGMECON method the model is as follow,: where is the RHS of the constrained objective functions, δ is a small number (usually ∈ 10 , 10 ).Also, in order to avoid any scaling problems, it is recommended to replace the (the surplus variable of the i-th objective function) in the second term of the objective function, by / , where is the range of the i-th objective function (as calculated from the payoff table).
In AUGMECON2, the improved version of AUGMECON, the objective function is slightly modified as follows (Mavrotas, 2013): This modification is performed in order to perform a kind of lexicographic optimization on the rest of the objective functions, if there is any alternative optima.For example, with this formulation the solver will find the optimal for and then it will try to optimize , then , and so on.With the previous formulation the sequence of optimizations of was indifferent, while now we force the sequential optimization of the constrained objective functions (in case of alternative optima).
In AUGMECON 2 like AUGMECON, for each objective function i=2…p the objective function range is calculated, then the range of the i-th objective function is divided into equal intervals, thus there would be total 1 grid points, that are used to vary parametrically the RHS of the i-th objective function.So the discretization step for this objective function is given as: The RHS of the corresponding constraint in the k-th iteration in the specific objective function will be as Eq. ( 64), where is the minimum obtained from the payoff table and k is the counter for the specific objective function: But in AUGMECON2, in each iteration, surplus variable that corresponds to the innermost objective function is checked.In this case it is the objective function with p=2.Then the bypass coefficient is calculated as: where int() is the function that returns the integer part of a real number.When the surplus variable is larger than the , it is implied that in the next iteration the same solution will be obtained and the only difference is the surplus variable which will have the value .This makes the iteration redundant and therefore it can be bypassed as no new Pareto optimal solution is generated.The bypass coefficient b actually indicate how many consecutive iterations can be bypassed.Thus, in this way, especially when we have many grid points, AUGMECON2 method can significantly improve the Solution time (Mavrotas, 2013).Accordingly, our proposed model in AUGMECON2, forms as follows: (66) where S is the feasible region of our proposed model.Also, the third objective function of our proposed model ( , is selected as the main objective function, and our first and the second objective functions are incorporated in its constraint part.Since our third objective function has fuzzy parameters, the expected value of it is used in the model Eq. ( 66).Other terms in the objective function of the above model, as described before, are related to solution robustness and model robustness of our proposed model.

Evaluation of the proposed robust model
First, we assume a natural disaster such as an earthquake happens in the center of three concentric circles, consisting of an inner, middle, and outer circles.Second, the AA nodes are randomly placed in the inner circle, as the RDC nodes, and the supplier nodes are randomly generated in the middle and outer circles, respectively.
Test problems are generated by randomly changing the values of some critical parameters such as the number of suppliers, the number of relief distribution centers, and the number of potentially affected areas.In addition we assume: Test problems are generated by randomly changing the number of suppliers, the number of RDCs and the number of potentially AAs.The number of suppliers is from 1 to 10, the number of RDCs is from 3 to 20, and the number of AAs is from 4 to 25.After locating the network nodes, parameters such as , , are estimated according to the Euclidean distance between these nodes and the earthquake epicenter, the disaster type, the disaster intensity, and urban fabrics.To generate the demand amounts, first the population of each AA is randomly generated.Then, water and food demands at each AA are estimated on the basis of the population, multiplied by the vulnerability probability ( ) of the AA.As the considered shelter has a capacity for accommodating three people, the estimated shelter demand at each AA is set to the number of affected people divided by three.Then, by the estimated demand amounts, we can generate the fuzzy trapezoidal demand numbers, as follows: Since there is usually sufficient time before disaster to supply an abundant amount of relief commodities, we estimate the pre-disaster total supply of each commodity to be about three times of its total estimated demand.The post-disaster supply amount is usually limited and uncertain, and we estimate it as follows: The above formula is used for generating the water and food amounts and one third of this amount is used as the shelter amount.Then, by the estimated post-disaster supply amounts, we can generate the fuzzy trapezoidal post-disaster supply numbers, as follows: To estimate the relief pre-disaster transportation costs, we multiply the inter-node Euclidean distances of our relief network to the unit transportation cost of each commodity.The post-disaster unit transportation costs are considered 1.2 times that of the pre-disaster.The fuzzy trapezoidal post-disaster transportation costs are defined as: 2 , 2 , 2 , 2 2 0.8 , 0.9 , 1 2 As proposed by Bozorgi-Amiri et al. (2013), holding cost of a unit commodity in response phase can be estimated according to its procurement cost, and the post-disaster procurement cost of each commodity can be considered nearly equal to that of pre-disaster.So the fuzzy trapezoidal post-disaster procurement costs are estimated as follows: Accordingly, to evaluate the effectiveness of the proposed robust model against its equivalent deterministic form, the following steps are performed: 1. Five different test problems are randomly generated (as described before).2. The problems are modeled based on both the robust and the deterministic forms.3. Models are solved by LP-metric solution method, using GAMS 23.9 running on a PC Pentium IV-5GHz with 4GB of RAM (DDR 3) under the Windows 7 environment.Moreover, for each problem, pre-disaster Strategic and tactical decisions, such as the location, size and retrofitting level of the RDCs, and the pre-positioned amount of commodities, are determined under both of its robust and deterministic models.4. Ten different random scenarios are generated for each problem.A scenario is a specific realization of the uncertain parameters of a problem.So, each problem will be modeled in a deterministic form under each of its related scenarios.5. Two surveys are conducted under each of the ten scenarios of a problem so that, each time, the pre disaster decisions adopted by solving one of the two models in step 3, are inserted as input parameters in the deterministic model made in step 4. 6. Deterministic models of step 4 are solved using LP-metric solution method and their result are used to evaluate the effectiveness of the robust approach.Furthermore, it is assumed that the importance of the first and the second objectives are respectively four and three times the importance of the third objective.Also, it is assumed that the relief organization is planed to meet at least 90% of AA's demand.In Fig. 1, the "robust OFV" and the "crisp OFV" are respectively representing the objective function values of the robust and the deterministic models generated in step 5.As Fig. 1 demonstrates, the robust OFVs remain almost constant in different scenarios over the five test problems, while the crisp OFVs fluctuate so much.Moreover, in some scenarios, the pre-disaster decisions adopted based on the deterministic models are unable to satisfy the second constraint of our proposed model (meeting at least 90% of each AA's demand), so the problem becomes infeasible and in Fig. 1 its crisp OFV column remains empty.As Fig. 2 presents, in pre-disaster planning, both the robust and the deterministic models of step 3, claimed a full demand response (Z1=0).However, in deterministic model the expected value of Z1 under 10 different disaster occurrence scenarios, shows that the deterministic model is incapable of a full demand response in post-disaster situations.While, as Fig. 2 shows, the average performance of robust models under those scenarios can prove their claim of a full demand response, and in this way, it proves the effectiveness of the robust model in after-disaster environment.
Also the results for the second objective function (Z2) show that the values of Z2 in robust models, unlike the deterministic models, were, in most scenarios, less than their pre-disaster prediction values.So, considering the high importance of the second objective function for the decision maker, this matter can improve the reliability of the robust model's pre-disaster decisions.Given the demand created by the recent earthquakes, it seems generation of a minimum 50% demand and a maximum 88% demand are suitable for its future earthquakes.Therefore, by considering a change interval, such as (0.8, 0.9, 1.05, 1.18) around the demand values of the 23 AAs, their trapezoidal fuzzy demand numbers are derived.Table 5 presents the information about the six selected suppliers in this case study.Table 6 contains the information of the relief commodities in pre-disaster phase.In the response phase, the unit procurement cost and the unit transportation cost are estimated to be respectively 1.1 and 1.2 times that of their predisaster phase counterpart costs.Also, to calculate the transportation costs all around the relief network, the distance between its nods is estimated using Google Maps.Finally the trapezoidal fuzzy numbers of the supply, demand and the after-disaster relief costs are obtained using proper change interval around their pre-disaster deterministic numbers.Moreover, due to rural characteristics and cold climate of this area, relief commodities priority order is assumed as: shelter, water and food respectively.In this section, the proposed problem is solved using the improved augmented -constraint method, by GAMS 23.9 running on a PC Pentium IV-5GHz with 4 GB of RAM (DDR 3) under the Windows 7 environment.To run the proposed model in GAMS, DICOPT solver with CPLEX as its MIP solver is used.In Fig. 6, Pareto optimal solution set for the proposed case study problem, Considering 30 separating points in the improved augmented -constraint algorithm and with 615 Pareto optimal solutions, is depicted as follows:

Conclusions
It has been the first time that a robust chance constrained programming model was provided for disaster relief logistics.We proposed a three-objective model that simultaneously considered the distributive justice in relief distribution, the risk of relief distribution, and the total logistics costs, which was then solved using LP-metric and Improved augmented -constraint (AUGMECON2) methods.This model not only considers the probability of facility (suppliers and RDCs) disruption during the disaster, but also as a novel work considers the disaster retrofitting of distribution centers (in preparedness phase).Also, relief commodities priority and affected areas (AAs) priority are considered in logistics relief distribution.In fact, disaster retrofitting and priority concepts in relief delivery are provided to improve the quality of relief services in the catastrophic disaster conditions.Moreover, the model is committed to meet a minimum percentage of each AA's demand level, which exerts the distribution standards of the relief organization throughout the network.Also, the uncertainty of relief supply in considered in a new fashion, by distinguishing between the suppliers' supply ability in preparedness and response phases, which helped us having a more realistic plan for relief logistics.Then the evaluation of some test problems showed that robust approach despite higher relief costs, can provide a good estimate of the crisis costs and can increase the reliability of the strategic decisions.This is while, deterministic approach could not even meet the minimum relief coverage standards of the relief organization in 38% of the considered scenarios.
Finally, a seismic region of northwestern Iran (East Azerbaijan) selected as a case study, due to its historical earthquakes, the large number of casualties and losses, and the rural characteristics of this region in recent earthquakes, and its relief logistics network modeled and proper strategic and operational decisions in the face of future earthquakes in this region reported.Which showed that the proposed model can be useful for decision makers in the field of natural disasters.
23)Aaccording toPishvaee et al. (2012), we can calculate the expected value of the fuzzy numbers of the above equation.In addition, constraints of the model are as follows, Constraint 1:

Fig. 2 .
Fig. 2. Evaluation of the first objective-function (z1) for both deterministic and robust models in five test problems

Fig. 5 .
Fig. 5. Potential AAs, potential RDCs and suppliers in the sixth seismic zonation of East Azerbaijan

Fig. 6 .
Fig. 6.Pareto optimal solutions of the East Azerbaijan relief logistics modeling

Table 1
Structure of facility location models based on the data type and number of levels and objectives Holding cost of one unit commodity c in RDCs, Required unit space for each unit of commodity c, The capacity of an opened RDC with size l (in cubic meters), Disaster risk index of RDC j based on its proximity to the center of the disaster, Importance factor of AA k based on its proximity to the center of the disaster, Importance factor of commodity c in the disaster, Distance between RDC j and AA k, An arc index for the path between RDC j and AA k, Big number in mathematical modeling, Integer numbers chosen by the decision maker , , Minimum percent of the demand of each AA that should be responded.Variables:Amount of commodity c transferred from supplier i to RDC j in response phase, Amount of commodity c transferred from supplier i and pre-positioned at RDC j in preparedness phase, Amount of c transferred from RDC j to AA k in response phase, notations used in our model are summarized as follows: Indices: Relief supply points, i Relief distribution centers or RDCs, j,e Affected areas or demand points or AAs, k Relief supplies, including food, water, and shelter, c RDC sizes, which include small, medium, and large, l RDC disaster retrofitting levels, including level 0 (not retrofitted), level 1, and level 2 (totally retrofitted).re Parameters: Building cost with size l (before the disaster), Retrofitting cost for an RDC with the retrofitting level of re, Unit transportation cost of one unit commodity c from supplier i to RDC j in preparedness phase, Unit transportation cost of one unit commodity c from supplier i to RDC j in response phase, Procurement cost of one unit commodity c from supplier i in preparedness phase, Procurement cost of one unit commodity c from supplier i in response phase, Unit transportation cost of one unit commodity c from RDC j to RDC e in response phase, Unit transportation cost of one unit commodity c from RDC j to AA k in response phase, Unit amount of commodity c supplied from supplier i in preparedness phase, Unit amount of commodity c supplied from supplier i in response phase, Unit amount of commodity c demanded at AA k, The percentage of stored amount of commodity c at RDC j with retrofitting level re, that remains usable in response phase , The percentage of commodity c at supplier i that remain usable in response phase, Shortage cost of one unit commodity c in response phase, Shortage amount of commodity c in AA k in response phase, Inventory amount of commodity c holding at RDC j in response phase, Amount of commodity c transferred from RDC j to RDC e, as a backup coverage, in response phase 1 if an RDC with size l and retrofitting level re is located at candidate point RDC j, 0 otherwise, 1 if RDC j sent commodity to AA k, 0 otherwise.

Table 2
Characteristics of the randomly generated test problems Evaluation of the LP-metric objective-function values for both deterministic and robust models in five test problems

Table 3
displays candidate RDCs and the approximation of for each of them.Also RDCs construction costs, according to their storage capacity, is presented in Table4:

Table 3
The candidate relief distribution centers for the case study

Table 4
RDC setup cost depending on its storage capacity

Table 5
The case study relief supplier points information

Table 6
Procuring price, transportation cost and volume occupied by each commodity unit in preparedness phase