A New Multi Objective Mathematical Model for Relief Distribution Location at Natural Disaster Response Phase

Every year, natural disasters such as earthquake, flood, hurricane and etc. impose immense financial and humane losses on governments owing to their unpredictable character and arise of emergency situations and consequently the reduction of the abilities due to serious damages to infrastructures, increases demand for logistic services and supplies. First, in this study the necessity of paying attention to locating procedures in emergency situations is pointed out and an outline for the studied case of disaster relief supply chain was discussed and the problem was validated at small scale. On the other hand, to solve this kind of problems involving three objective functions and complicated time calculation, meta-heuristic methods which yield almost optimum solutions in less time are applied. The EC method and NSGA II algorithm are among the evolutionary multi-objective optimization algorithms applied in this case. In this study the aforementioned algorithm is used for solving problems at large scale.


1-Introduction
Natural disasters cause, millions of casualties and leads to reduction of capabilities and financial losses across the world every year. As a developing country Iran is one of the geographical regions most prone to unpredictable emergency incidence and it is known as one of the ten most susceptible countries to natural disasters which almost 90 percent of its population are exposed to natural disasters. Considering the impact of natural disasters on the

International journal of innovation in Engineering
journal homepage: www.ijie.ir community's health, wellbeing, and welfare, providing healthcare services is the primary factor of survival and mortality reduction and individuals' welfare in the next phase after the occurrence of such incidents. Needless to say when encountering crisis and disasters, the initial action that will be taken is rescuing and mitigating the effects of the incident which requires quick response considering the limited time. Quick response which constitutes substantial part of crisis management, involves identification, assessment, decision making and provisional emergency measures, which whole part of this procedure takes place within a very short time sometime a few hours, therefor one of the actions that will be taken in order to manage the crisis is to consider relief measures aftermath of the crisis. The humanitarian environment of human rights organizations are built on three principles of humanity, impartiality and fairness (justice). In other words, they help any person who is in need of help, at any place regardless their beliefs even though in contradiction to theirs and do not privilege any group of final recipients over other beneficiaries (Webber, 1999). Humanity means that all the suffering and pain should be eliminated wherever it is. Fairness indicates help with no discrimination prioritized based on emergence of the needs. From physical point of view, the humanitarian environment indicates a peaceful area in which civilians and non-militant individuals and non-armed relief workers are able to freely move and operate. Often it is difficult to establish humanitarian principles in the complicated environments, especially in the situation of war and military conflict (Webber, 1999). Humanitarian logistics usually operates under a high level of uncertainty regarding demand and required supplies. Humanitarian workers often deal with many beneficiaries including aid workers, intermediaries, governments, and military beneficiaries and final recipient in a logistic humanitarian environment. The lack of coordination among humanitarian organizations in the crisis situation is common, as such there can be hundreds of humanitarian organizations with different political instructions, ideologies and religious beliefs aiming to carry out relief operations. The biggest challenge in these situations is the levelization of these organizations regardless of their beliefs and authorities (Farahani Hekmat-Far, 2009). The supply chain management concept is based on this concept that supply chains can outperform individual business units in the competition. A supply chain is consisting of a range of beneficiaries including suppliers, manufacturers, distributors, retailers and customers. The supply chain's structure varies for different companies even within similar industries. If the commercial supply chain is defined as a network streaming products, information and financial resources from the source to final customers, the relief or humanitarian supply chain therefore can be defined similarly as a network for managing the supplies, information and financial resources flow from relief sources to people affected by the crises. In the contrary to commercial supply chains which are owned by stockholders, the humanitarian supply chains are not centered around a targeted profit and they rely on help of volunteer forces and aid workers. Various models have been developed for establishing efficient supply chains based on cost minimization (which is equivalent of profit maximization) which can be used for establishing this kind of supply chains in their original forms or after some modification. An effective humanitarian supply chain management is capable of responding to many actions often on a global scale, as fast as possible and within a short timeframe.
Thereby humanitarian supply chains must be multifold, universal, dynamic, and provisional in nature. Relief logistics has many actors, including a large number of private sectors, donors, and host governments and the military, governments of neighboring countries and NGOs, all of which have different political, ideological plans and religious attitudes, different levels of media attention and donors (Figure 1). When a crisis or a disaster occurs, depend on the location, type and intensity of the disaster different items and supplies are required. It is necessary to categorize the relief supplies due to their large quantity and variety which these categorizations will be very useful in identifying and evaluating required items and as well as designing the distribution system (Bozorgi Amiri et al., 2013). In one of these categorizations, relief items are segmented into two main categories of consumable and non-consumable ( Figure 2).
As demonstrated in Figure 2, the consumable relief items are distributed to the people affected by the disaster for several times, while non-expendable relief items require one-time distribution.

Relief Government
Other NGOs

Logistic Suppliers Donors Military
Disaster Relief Logistics is one of the main activities in the subject of disaster management. The Relief Logistics means to distribute the right materials to people in the right places at the right time and to the right amount (Louveaux, 1993). The approaches applied include two steps: in the first step (preparedness phase), the location of relief distribution centers, the amount of inventory of relief supplies for storage purposes and for appropriate suppliers are determined. the next step (response phase) include decision making regarding the appropriate suppliers in the crisis occasions, the amount of relief supplies distributed from supply points to relief centers and from these centers to the areas affected by the disaster. Aid distribution and required services in a timely manner is one the important subjects in relief logistics (Bashiri et al., 2009). The summarized relief logistics operation is demonstrated in figure 3 (Bremberg & Joyale, 2001).

Literature Review
Natural disasters and man-made disasters such as floods, earthquakes, hurricanes, droughts, wars and etc. affect people in different parts of the world every year. These natural crises (such as earthquakes, floods, hurricanes) owing to their unpredictable character cause immense and sometimes irreparable damage to nations. In fact, more than 599 disasters occur on our planet every year, kills nearly 75 thousand people and affect more than 200 million people (Bremberg & Joyale, 2001). For example, nearly 1,836 people lost their lives in the Katrina hurricane in 2005, and at least 86,000 people were killed in the Kashmir earthquake in the same year (Rodriguez et al., 2004;Levine and Arisil, 2004). Also, Haiti earthquake in 2012 affected nearly three million people. These natural disasters often cause huge financial and humane damages at time of occurrence. The top five disasters with the highest human casualties during a period from 1980 to 2012 is available in table 1 (Camacho Vallejo et al., 2014). Unfortunately, the natural and man-made disasters will increase to five times for the next 59 years, hence the need for relief services will rise as well. Thereby, it can turn the humanitarian aid delivery into a global industry. Most of the time data in disaster management are uncertain, and quick decisions are made based on just a fraction of the required information. Not only technology but also social, cultural and psychological factors play an important role in decision making (kumawala, 1972). Hence upon, this study proposes a Robust Possibilistic Programming, in crisis situations under uncertainty.   Altay (2006) proposed a model for crisis management centers location for relief supplies management. Jia et al. (2007) reviewed facility locating models for a large-scale emergency situation and categorized it into three groups: coverage models, median location and center location models. Yi and Kumar (2007) proposed a meta-heuristic algorithm (Ant colony optimization algorithm) to solve the disaster relief logistics problem. Balcik and Beamen (2008), proposed a model for determining the number and locations of distribution centers in rescue operations. They formulated the location problem of another type of maximum coverage problem for a set of probable scenarios. The proposed objective function maximizes the total expected demand covered by the distribution center. They also determined the amount of relief supplies required to be stored at each distribution center to meet the demand. Although they did not consider locating as part of the supply chain, their study is one of the first facilities locating problem in rescue operations that must be solved. Braldi and Bruni (2009) proposed a possiblistic model for determining the optimal location of emergency facilities under uncertainty. Table (3) demonstrates a summary of the proposed models by different researchers.

Research Methodology
In this section, we propose a multi-objective mathematical model for the locating-distributing problem of relief items in the disaster response phase. • Yks If, 1 in the candidate's place, k is facilitated by the surface s and zero is otherwise.
• Yjklm Equivalent to the amount of cargo l of the product, using the transport m from the storage center j to the k region

• Zj
If the relay distribution center j is selected to be opened, it is equal to 1 and otherwise 0 • When a product l is required to receive a shipment of m from a donor country or an international organization i to a storage facility j.

•
The time needed to distribute a ll product by using m transport from the storage center j to the damaged area k • Ddjs The share of demand for the surface s of the demand point d, which is met by the distribution center j • Ddks A share of the demand level s of the demand point d that is satisfied by the affected center k.
• Djks A share of s level demand, which is referred from the distribution center j to k (ie, the level s of the demanded amount of the injured) • Yks If, at the candidate location, k is the center of the distribution of the surface s, then zero is set to zero otherwise

• Fj
It is equal to one when the distribution center is established at j, otherwise zero Gets • Is equal to one that is assigned the type t at the point of demand d to the damaged point k , Otherwise it becomes zero.
The objective function is described as follow: Objective 1 minimizes the total cost of supplies shortage.
Objective 2 is the total cost of delivering demanded supplies and the cost of constructing the distribution center, which includes four group of costs as follow: the first group of delivering demanded supplies between customers and the existing distribution center, second group is delivering demanded supplies between customers and the recently established distribution center, the third grope is delivering demanded supplies between distribution centers and the areas affected by the disaster. The fourth group consists of fixed expenses.
Objective 3 minimizes the total response time which is defined as the required time to deliver the aid.

The procedure of presenting a model for demand uncertainty
According the presented model in this section, different fuzzy approaches to planning are considered.

BPCCP model
For modeling the objective function, the operator's expected value and the appropriate value that will fulfill the possible constraints with vague parameters were used. Also, we considered the trapezoidal probability distribution with four outstanding points as = ( (1)

RPP model
First, we assume that the f, c, and d vectors are ambiguous parameters, and in this phase the coefficient matrix N are not vague. Therefore, the BPCCP model is defined as the RPP model. This model aim to establish balance between the sections of the objective function: 1. average performance 2. optimality robustness

Feasibility robustness
In this model same as the BPCCP model, the first term in the objective function indicates the expected value for z and the second term indicates the difference between two final possible value for z which are defined as follow: = 4 + 4 (21) And parameter indicates importance (weight) of this term against other two terms of the objective function. Thus, the second term ensures that the gap between the maximum and minimum expected optimal value Z of the objective function will be obtained. In fact, this term controls the optimality robustness of the solution vector. The third term determines the confidence level of any possible constraint, and the penalty for any possible violation of any constraint, including ambiguous parameters, and shows the difference between the worst possible value of the uncertain parameter and its value in the chance constraints; in fact, it controls the feasibility robustness in the solution vector. It should also be noted that it is not only a theoretical or meaningless parameter, but it can be determined based on application in the right context.
: , , , , , , ∈ As you can see in the model, when technical coefficients are assumed uncertain, the linearity of the constraints and the possible objective function will turn to nonlinear, although in some cases a nonlinear term can be turned into linear with defining a new variable and adding few constraints to the model. To avoid the complexity of the nonlinear model, the auxiliary variable and the nonlinear model will be turned into an equivalent linear model.
Thereby, the above nonlinear model will turn into the following linear model:

Results and Discussion
In this section, a solution for the problem is discussed and furthermore, analysis and interpretation of the results are presented. The NSGA II meta-heuristic algorithm was used to solve the problems, and the computer programming was conducted with MATLAB R2010a and was carried out using a portable computer with CPU: 2.30 GHz and 4 GB of RAM. Also, to generate different values for the two objective functions that have been defined in this problem a fuzzy numerical simulation algorithm was applied, which is described as follow. Furthermore, a neural network system with a fuzzy inference engine is designed for the problem.

Problem Solving Method
Given the multi objective nature of this mathematical model, a NSGAII Meta heuristic algorithm is used and the results are presented. Also, to achieve an efficient objective function, we designed a fuzzy neural network to find the Pareto solution, which will be fully discussed in this section. It is worth noting that a fuzzy simulation algorithm was used to generate different values for the objective function problems.

The NSGAII Algorithm
The NSGAII algorithm which uses non-dominant sorting is one of the most well-known and widely used optimizing algorithms for multi-objective optimizing. This algorithm was first proposed by Deb in 2002. Alongside all the application that is provided by the NSGAII has, it can be considered as the formation pattern of many multi-objective optimization algorithms. This algorithm and its unique approach to multi-objective optimization problems have been frequently used to develop new multi-objective optimization algorithms.

Comparative Metric for Evaluating Multi-objective Algorithms
There are two main metric groups to evaluate the performance of multi-objective meta-heuristic algorithms, including convergence metrics and dispersion metrics.
The first metrics group is number of Pareto solution and Mean of Ideal Distance, the second group metrics are spacing and Diversification and also the time metrics. According to work of Bashiri, 2009 andAfshar, 2011, four metrics was applied to compare these algorithms in this study, which are explained in the future.

The Mean of Ideal Distance Metric
This metric is used for calculating the average distance from Pareto solutions to the origin of coordinates. According to this equation, it is obvious that the lower this metric is, the greater the efficiency of the algorithm will be. Since, one of the goals of multi-objective Pareto-based approach is to have far closer solution to the origin boundaries, thus this metric calculates the distance of the fronts from the best value for the population. The n is the numbers of Pareto solutions and the mean distance from ideal points is obtained as the value for this metric. Because the type of objective function is minimization with equal weight for each objective, we selected the origin point as a best value.
In this thesis, considering the fact that the both objective functions are minimization functions, thus the ideal point is  (31)

Algorithm running time metric (CPU T)
In the large scale problems one of the important indexes is the solution running time, therefore algorithm running time is considered as one of the quality assessment indexes.

Spacing metric (SM)
This metric represent the uniformity of Pareto solution distribution in the solution settings. This metric is calculated as follow: The lower the SM index, the higher performance of the algorithm.

Diversification metric (DM)
The diversification metric measures the diversity range of the Pareto solutions of an algorithm and it is calculated according the following equation. The higher the DM index, the higher performance of the algorithm.

Experiment design to set algorithm parameters
As it was discussed in the previous sections, the NSGAII algorithm involves many factors and parameters, affecting the final results and performance of the algorithm, therefore the algorithm performance is vastly affected by the appropriate combination of these factors. One of the objectives of experimental design is to observe and identify the shift in the output variables corresponding to intended shifts made to the input variables.
There are several methods to design an experiment. One of the first methods presented in this field is the factorial There are two methods for analyzing the experiment: 1-Standardized method of ANOVA (Analysis of variance)

2-use of signal-to-noise ratio (S/N)
The Taguchi method classifies the objective functions into three categories: the lower the better, the higher the better, and the best nominal. Since our objective function is intended to minimize the prediction error, it is considered to be in the first categories. In this research, the proper values for the algorithm parameters are determined during the sequential steps with the help of Minitab software.  Now, in order to determine the optimal parameters, the criteria are converted to a solution using one of the multiple attribute decision making methods, called Simple Additive Weighting (SAW). Accordingly, the instance problems were generated and solved for all of the 22 levels, and the results were imported to the Minitab software for the Taguchi experimental design. The settings of parameters resulted from solving the instance problems and implementation of the aforementioned steps are depicted in the figure 6. The optimum combination for the parameters of the NSGAII algorithm can be set by examining Figure 7. Table 5 shows these values. According to table 5, values including 0.9 for crossover probability, 0.2 for mutation and 150 for initial population were set. Also, according to the initial experimental data, the number of generations required to improve the algorithm is set as 50.

Problem Solving using the GAMS software
The operational research software, GAMS is a high-level programming language used to solve mathematical programming models. This software is a very high-speed problem solver used for large scale models. In fact, GAMS can be considered as the best software for solving the most sophisticated and large scale optimization problems. The GAMS software easily connects to databases. Also, if problem doesn't involve high volume of input data, then the data can be entered into the software directly. For modeling, all variables and constraints are set once. Since the software easily iterate restrictions and there is no need to insert all the restrictions. The GAMZ software output is easy to use. The outputs can be specified in the software itself, or it can be for example, imported to the Excel after connecting to the databases. In this thesis 22 problem instances are solved using GAMS software, which the results of solving the problems are shown in the table (6).

Random instance generation
In order to verify the proposed model and the efficiency of the proposed algorithm, five problems are designed in different dimensions, which their related information are demonstrated in table (7). Also the specified bounds for the problem parameters are demonstrated in table (8). In table (7), the first column shows the problem numeral, the second column shows the number of supply points (first level), the third column shows the number of distribution centers, the fourth column shows the affected areas (second level), the fifth column represents the number of relief supplies, the sixth column shows the amount of demand, and the seventh column shows the quantity of demand, column eighth shows the number of transportation vehicles, an Eighth d the ninth column represents the type of equipment.  (8), the required parameters of this problem were randomly generated using uniform distribution at proper time intervals. After the random generation of the instances in various dimensions, the problems were implemented in the GAMS software, and these problems were ran in different sizes with proposed values, and the result are demonstrated in the following tables.

Examining the proposed NSGAII algorithm efficiency and solving the problem instances
In the proposed ε-constraint Method, the first objective function is considered as the primary objective function and the second one as the secondary objective function; next, 10 breakpoints were specified for the second objective function, and a total of 10 Pareto points for each problem is generated. In this section, the five-dimensional problems designed in the previous sections, notwithstanding both objective functions, are solved using the ε-constraint Method and the proposed NSGA II algorithm, and the resulted Pareto solutions are reported. In order to explain the ε-constraint Method, we consider the No. 3 problem. First of all, in order to set the Pareto points based on each objective function the problems are solved separately. The results are reported in Table ( (7) to (9), the resulted Pareto points using the two algorithms for problem No.3 are demonstrated as an example in small scale. It is noteworthy that these points and objective functions are shown correspondingly.       (12) to (16) and graphs in figures 10 to 13 it is evident that in small scale (problem 1 to 3), to the extent that exact solving the problems in 3200 seconds timeframe is possible, the results for NSGA II algorithm and the ε-constraint Method are very close. Accordingly, the metaheuristic algorithm used in this study results in solutions with close proximity to the exact solutions, and thereby it can be considered as a proper tool for problem solving in the cases that exact solutions are ineffective. For example, starting from Problem No. 4, the ε-constraint Method is not capable to solve problems within the timeframe exactly.
As such, the NSGA II algorithm is used for solving problems on the same aforementioned scale, on account of the adequate performance of this algorithm. Figure 10: The values of mean distance from the ideal solution metric for the two algorithms. Figure 11: The values of spacing metric for the two algorithms. Figure 12: The values of diversity metric for the two algorithms.
The graphs of the spacing metric and the mean distance from the ideal solution, distance and diversity metrics, leads us to same conclusion that the NSGA II algorithm has an adequate performance for all of the metrics and close proximity to the exact method, which accordingly its final performance will be as close. The final comparison of the two methods is shown in figure 14 which confirms this conclusion.   (17), as the scale of the problem grows the required time for exact solution of the problem increased dramatically to the extent which after problem instance No.3, the ε-constraint Method is unable to keep up with the timeframe set finds exact solution of the problem. Instead, the meta-heuristic algorithm is able to solve the problem in much shorter time. Therefore, the NSGA II algorithm has a justifiable performance during the proper timeframe set for solving the problems. It can also be observed in figure 14. Figure 14: Timeframe for both EC and NSGA II solving methods (seconds).
Ultimately, given the adequate performance of the NSGA II algorithm, large-scale problems in this research are solved using this algorithm and, as an example of large-scale problems, the optimal Pareto-frontiers resulted using the NSGA II algorithm for problem No.5 are presented in figures 15 to 17. Henceforth, the definition of the NSGAII algorithm with non-dominant sorting, an evolutionary multi-objective optimization algorithm, was discussed. In NSGAII algorithms and NSGAII based algorithms, each chromosome represents a point in the solution space and it is a feasible solution to the problem. To demonstrate a feasible solution in the proposed algorithm, a string of one and zero with the length of I*K is used, which I represent the same number of activities and K is the number of modes for each corresponding activity. And in case of selecting the mode k for activity i, the bit associated to that chromosome is zero and, otherwise it is one. Also the fitness function for the objective functions was defined. Considering 3 objective functions simultaneously for the problem, after cost calculation in the string was performed for once, then time and cost and distance constraints were imported into the three objective NSGAII algorithm, the Pareto solution are shown in the tables, which objectives contradiction is evident. In another method, three objectives counterbalance were solved using NSGAII algorithm and Pareto solutions associated to each solution are presented. And for each mode, the lowest cost and time is calculated.

Conclusions and Future research
In this study, first the necessity of paying attention to locating procedures in emergency situations is pointed out and an outline for the studied case of disaster relief supply chain was discussed and the problem was validated using GAMS and Matlab software at small scale. On the other hand, to solve this kind of problems involving three objective functions and complicated time calculation, meta-heuristic methods which yield almost optimum solutions in less time are applied. The EC method and NSGA II algorithm are among the evolutionary multi-objective optimization algorithms applied in this case. In this study the aforementioned algorithm is used for solving problems at large scale.
One of the drawbacks with this model, is the required complicated calculations, and among its advantages is that it involves three objectives, including the cost of relief supplies shortages, costs of the supply's transportation, and the costs of constructing new distribution centers and finally the response time to access to people affected by the disaster which is innovative in compare to other models. The directions that can be considered for the future research are as follow: