Multi-level facility location-allocation problem for post-disaster humanitarian relief distribution: A case study

Sets

V

the set of nodes (v, v′∈V)

J

the set of candidate locations (j, j^'∈J, J⊆V)

Scalars and parameters

M

A large positive number

d_jj′

Euclidean distance between facility j and facility j′

d_jv

Euclidean distance between facility j and node v

λ_υ

the demand located on edge (v, v′)

E

the endurance of each drone

M_t

the maximum allowed length of the route from demand point to its assigned launching station

k

the maximum number of relief centers to be established

N

the maximum number of recharge stations to be established

Variables

L_j=1

if a launching station is established in location j; 0, otherwise

R_j=1

if a refueling station is established in location j; 0, otherwise

y_j=1

if a facility is established in location j; 0, otherwise

x_jv=1

if node v is assigned to facility j; 0, otherwise

t_jv

the length of the shortest path from facility j to demand node v

• A multi-level FLP is used to simultaneously account for both relief centers and recharge stations.

• The model aims to optimize waiting times, total travel distance, survival rate and topology of the humanitarian relief system.

• The number of required resources including centers and drones is minimized without risking the efficiency of the system.

• A hybrid genetic algorithm is proposed for solving the model.

• A case study is investigated to evaluate the performance of the proposed method.

Millions of disasters occur each year on the planet Earth. Only during 15 years between 2000 and 2015, there have been 800,000 lost lives just because of earthquakes (Lamb and Jones, 2012). In years 1999–2008, 1.2m lives were lost due to 7,000 disasters worldwide (Overstreet et al., 2011). Each year disasters disturb the lives of 250m individuals (United Nations, 2007). Recent incidents such as the tsunami in the Indian Ocean in 2004 and the 2003 storm in the Philippines have further highlighted the significance of humanitarian relief (Garner and Harrison, 2006; Carlson et al., 2016). The trend of disaster and rate of their occurrences have been positive in recent years and vulnerability is increasing due to increased population and ill-defined disaster response infrastructure (Apte et al., 2016; Soneye, 2014). The number of disasters is forecasted to be fivefold within 15 years (Thomas and Kopczak, 2005) and this has made governments have more emphasis on disaster management in their agenda (Kovács and Spens, 2007).

Disaster management includes five phases: prevention, preparedness, response/relief, rehabilitation and reconstruction. If planned properly and deployed promptly, disaster response can significantly decrease the fatality rate (Khalil et al., 2014). Thus, preparedness and response phases of disaster management are most discussed, while other phases specifically recovery phase has not been addressed as much (Leiras et al., 2014). Furthermore, mitigation, preparedness, response and reconstruction are counted the most vital phases of disaster management (Cozzolino, 2012; Jahre et al., 2016). Relief distribution is a key factor in disaster response (Ozen and Krishnamurthy, 2018). Modern technologies should be utilized to distribute relief among survivors immediately after a disaster to increase the survival rate (Eguchi, 2013). Medicines, food, water and shelter are the most important elements of emergency relief (Safeer et al., 2014). Most important tasks of post-disaster relief are the proper allocation of relief crews and resources, locating the relief centers, and generation of transportation routes (Chang et al., 2014). Consequently, 80 percent of all humanitarian relief operations are built upon logistics without which the relief distribution would be corrupted (Chomilier et al., 2003). There are also many studies related to strategies and available technologies in disaster management and emergency response system and their discussion is out of the scope of this study (Delmonteil and Rancourt, 2017; Serrato-Garcia et al., 2016; Hong et al., 2015; Jahre, 2017; Makepeace et al., 2017; Vybornova and Gala, 2016; Gavidia, 2017).

Facility location problem (FLP) has a great importance in the strategic planning of organizations’ logistics. The wide spectrum of FLP applications has put it under the spotlight in many studies in the context of supply chain management (SCM), healthcare management, emergency logistics and so forth. Although SCM was introduced as an independent field, OR and FLP found their way into its core (Melo et al., 2009). SCM is the process of planning, implementing and controlling the operations of the supply chain in an efficient way. Part of the planning processes in SCM aims at finding the best possible supply chain configuration and that is where FLP plays its role in SCM. The interested reader is referred to Melo et al. (2009) for a complete review of the relation between SCM and FLP. In the case of humanitarian relief, the design of supply chain significantly affects its performance (Hong et al., 2015). It has been indicated that improper selection of the facility locations results in poor accessibility and functioning of the facilities (Rahman and Smith, 2000). Facility location received importance due to its direct impact on operational cost and timeliness of response to the demand (Haghani, 1996). Thus, many studies have investigated the proper location of different units and elements of a humanitarian relief system (Maharjan and Hanaoka, 2018; Timperio et al., 2017). In his study, Tayal and Singh (2017) reviewed the application of FLP for disaster relief operations. Generally, FLP aims to determine the best location of new facilities among a set of candidate locations. Facility location models for profit-making bodies is primarily to minimize cost or maximize revenue, whereas facility location models for public and emergency services focus mainly on accessibility and response time. FLPs consider a facility opening cost corresponding to each location. In k-median FLP, there is no cost associated with facility opening and just total service costs or aggregate travel distance is minimized. The k-center FLP is defined as the minimization of the maximum travel distance. In the covering problems, all the demand is supposed to be satisfied and there is no lost demand. Set covering-based facility location models select facilities among potential locations such that all demand sources are covered with a minimum number of facilities. Therefore, in disaster relief operations, this means that each potential demand location must be within the targeted response time of a facility in the relief network. However, it was observed that a small portion of the demand located in remote areas has a huge effect on the final solutions. Thus, Charikar et al. (2001) addressed the lost demand, called outliers, in the FLPs. Another categorization of FLP is based on the capacity of the facilities (Charikar et al., 2001). In un-capacitated FLP, it is presumed that there is no limit on the services provided by facilities and consequently no lost demand, while in capacitated FLP each facility has an upper limit for its services. As described by Laporte et al. (2015) in his book Location Science, fixed-charge FLP deals with a finite number of customers and a finite number of candidate locations and the costs of the facilities are considered in the model. Furthermore, there are generally two types of decisions made in these problems: determining the location of service providers and allocation of the demand. Fixed-charge facility location also plays a critical role in many other areas like SCM, distributed systems, humanitarian relief, emergency systems, location-routing problems or freight transportation. Traditionally, the facility layout design problem is linked with optimizing the performance of the plant by minimizing the total flow of the products. Therefore, the facility layout problem is dependent on the product flows and its demand. During a disaster, the demand may shoot up and the layout which was designed initially without considering the disaster may not be optimal now to handle the post-relief operations. Therefore, disaster criteria must be considered to design the facility layout, once the facility in a supply chain network is selected, to optimize post-disaster relief operations. Some scholars have used stochastic distributions to estimate the fluctuations in the demand during the post-disaster period (Tayal and Singh, 2017). FLPs are generally NP-hard on a general graph and thus exact methods cannot be utilized for solving them, especially when dealing with large-size problems (Laporte et al., 2015; Megiddo et al., 1983; Kariv and Hakimi, 1979; Contreras et al., 2011).

Drones have been widely investigated as delivery vehicles and due to its speed and capability of autonomous operation, it can be considered for relief distribution. Previously the use of drones has been investigated for humanitarian relief (Tatham, Ball, Wu and Diplas, 2017; Comes et al., 2018; Tatham, Neal and Wu, 2017), rescue operations (Xiang et al., 2016), delivery (Shavarani et al., 2017) and healthcare (Scott and Scott, 2017; Tatham, Stadler, Murray and Shaban, 2017). Golabi et al. (2017) considered the application of drones for the relief distribution and proposed a model to locate the best locations for the launch stations/relief centers. In their work, the delivery system was composed of aerial and ground transportation wherein the drones were used for the areas which were made inaccessible due to road collapses. On another hand, Hong et al. (2017) proposed a model for finding the best locations for recharge stations, assuming that there is only one launch station established. Only a few studies have scrutinized aerial delivery systems, and the authors have tried to cover and review most of the existing ones. Furthermore, there is no study about the logistics of the system in the current literature where both launch and refuel stations are considered, specifically in the humanitarian logistics and in the post-disaster environment. The current study conveys such a research by proposing a multi-level FLP and solution methods.

In this study, the drone delivery system consists of both relief centers and recharge stations which can be used by drones in order to reach to the points that geographically are not located within their endurance. To do this, a multi-level fixed-charge FLP is proposed to model the system and find the optimum number and location of both relief centers and refueling stations. Due to the NP-hardness of the problem, a genetic algorithm (GA) is devised to solve the model. The proposed algorithm is applied to the case of Tehran for finding the best center locations as well as refueling stations to be used for post-disaster relief distribution.

The structure of this paper is as follows: the problem is defined and formulated in the next section. In Section 3, a GA is devised to solve the proposed model. A case study of Tehran is presented in Section 4. Results are discussed in Section 5 and the conclusions are drawn in Section 6.

The set of candidate locations is assumed to be known. Among the set of candidate locations, a number of locations should be selected for establishing relief centers and refuel stations in a manner that the aggregate travel distance is minimized. It is assumed that the demand occurs according to Poisson distribution and is satisfied by its closest facility. Drones are identical and their speed is fixed, which implies that minimization of the total travel distance is equivalent of minimization of total service time. Each drone has an endurance that indicates the maximum distance it can travel without a need to refuel. If a demand were within the endurance of the drone, the drone would go directly to the demand point, deliver the relief package and turn back. Otherwise, to reach the demand point, the drone should visit one or more refuel stations. Thus, the length of the shortest paths between the stations and the demand points are not parameters but variables that depend on the opened refuel stations at each instance of the problem. These variables are solved for within the model by the constraints of the shortest path algorithm. The words refuel and recharge stations are used interchangeably throughout the text. The same goes for the launch station and relief center. The notations provided as nomenclature are used in model formulation.

Using the described assumptions and definitions, the problem can be modeled as follows:

Subject to:

The objective function minimizes the total travel distance to meet the demand on each point. Constraint (2) ensures that either a refuel station or a center be established in one candidate location. Constraint (3) guarantees that each node is assigned to one center for receiving services and constraint (4) prohibits the assignment of any node to closed centers. The maximum number of relief centers and recharge stations is limited by constraints (5) and (6). Constraints (7)–(9) provide the shortest paths between nodes and centers in a way that the endurance capacity of the drones is not violated. Constraint (10) provides an upper limit for the shortest paths. Constraint (11) makes sure that each node is assigned to its corresponding closest center. Constraints (12) and (13) define the binary and positive variables.

Due to NP-harness of FLPs (Owen and Daskin, 1998), metaheuristic algorithms are widely used to solve them. GAs are among the most frequently used algorithms in the literature (Maghfiroh and Hanaoka, 2018; Shavarani and Vizvari, 2018). GA provides competitive solutions with less CPU time compared with other algorithms used in this area such as particle swarm, ant colony and simulated annealing (Arostegui et al., 2006; Youssef et al., 1998, 2001). Holland (1992) developed GA based on the evolution of superior individuals through generations. This algorithm starts with a population of random solutions. Through each iteration, some members of the population are randomly selected as parents and by mutation and crossover operators, the selected parents are mated and offspring solutions are produced and added to the pool. The solutions are then sorted based on their fitness value. The next generation evolves by selecting the superior solutions.

As indicated in Figure 1(a), each solution is represented by a chromosome wherein the ith gene indicates the status of the ith candidate location. There are three different values for each gene: 1, 2 and 3 which, respectively, indicate relief center, refuel station and closed status. The initial population is consisted of randomly generated solutions. In order to generate offspring/children, parents are selected by roulette wheel technique that gives higher priority to superior individuals. The superiority of the chromosomes is determined based on their fitness value as defined in Equation (1). The less the fitness value, the better is the solution. Parents undergo the crossover and mutation operations to generate children. To perform a crossover operation, two parents are selected, a crossover point is generated randomly and the second part of the chromosomes is interchanged, giving birth to two new children, each inheriting some characteristics of their parents. For each mutation operation, only one parent is selected and some randomly selected genes would change within the chromosome of the selected parent and an offspring is generated. The chromosome structure, crossover and mutation operators and pseudocode of the algorithm are illustrated in Figure 1. To improve the quality of the solutions and to speed up the convergence rate toward better solutions, a local search is also applied to some randomly selected solutions. In each selected solution, some genes are changed to a lower level, and the quality of fitness function is assessed, in the case of improvement the old solution is replaced by the new one. In the next step, population, which is now the union of original population and generated children, is sorted according to the fitness value and the best solutions are selected to enter next iteration of the algorithm and the rest of the solutions are discarded. This is compatible with the nature in the regard that the superior individuals have a better chance of survival. The algorithm goes on for a predetermined number of iterations. The number of solutions in the initial solution, maximum number of iterations, the rate of offspring produced in each iteration by crossover and mutation operators and the rate of local search are parameters that are determined by the decision maker. In this study, population size, maximum number of iterations, the rate of crossover, the rate of mutation and the rate of local search have been set at 1,500, 50, 0.8, 0.5 and 0.2. Taguchi method is used to tune the algorithm and find the best values for each of these parameters. The details of Taguchi methods can be found in many studies such as Maghfiroh and Hanaoka (2018), Shavarani and Vizvari (2018) and Nejad et al. (2018).

Tehran, the capital of Iran, is selected as the case study of this study. This city is under a severe risk of a major earthquake in the near future (Yaghmaei-Sabegh and Lam, 2010); therefore, extreme precautions must be taken into consideration for dealing with this issue. Thus, this case study investigates the best locations of relief centers and refuel stations to be used in such an event. As illustrated in Figure 2, the open street map of Tehran with 124,171 streets and 92,471 nodes was acquired from ArcGIS 10. Among existing nodes, 182 potential locations have been selected. The proposed algorithm is coded and run in MATLAB (R2015a). The information about the existing faults and intensity of the corresponding earthquakes and destruction factor (β_i) is acquired from the reports provided by Japan International Cooperation Agency (JICA, 2000). The demand at each node is bases on destruction probability obtained as p_i=1−β_iPGA_i (Salman and Yücel, 2015), PGA_i is the peak ground acceleration at node i defined by Panoussis (1974) as PGA_i=α(e^0.8μ/(r_i+40)²) where r_i is the distant of node i to the activated fault and μ is the magnitude of the earthquake. The demand is calculated by the assumption of ray fault activation with the magnitude of 6.5 and the total demand is assumed 82,914. The maximum number of centers and recharge stations is assumed 5 and 10, respectively. The endurance of the drones is assumed 16 km. Table I presents the results for different values of K and N. The results for the second instance are depicted in Figure 2, where the maximum number of centers and recharge stations are equal to 7 and 10, respectively. The total travel distance of this solution is equal to 393,682 km which translates to 4.7km per demand, considering the average speed of 80 km/h for a typical drone, the average of each mission could be estimated at 3.5min.

Disasters are inevitable and consequently each year millions of lives are lost to them. However, the total number of casualties could be significantly decreased by proper disaster management. To enhance the efficiency of relief actions, disaster managers try to utilize the state-of-the-art technologies. In recent years, use of drones has been under the spotlight because of the provided speed and accessibility to remote and disconnected areas; both of which are vital for saving the lives of the victims. To use the drones efficiently, the required infrastructure should be managed properly to make the most out the available resources. This study proposes a tool to determine the topology of the stations, both refuel and launching stations, in a way that total traveling distance is minimized. Since the speed of the drones is fixed, this objective corresponds to the minimization of response time and increasing the chance for maximizing the survival rate. The case study of Tehran shows that the model can act as a decision support tool for determining the optimum number of stations and their corresponding locations and simulation of the consequences of each combination of parameters such as the maximum number of refuel stations and launching stations. The case study is investigated for different values of maximum number of stations and the results are satisfying since the waiting time for receiving services is limited to several minutes. Disaster response plays an important role in disaster management and prompt response would significantly decrease the fatality rate. Findings of this study address the logistics of an aerial delivery system as a respond to the universal attention toward using technological advancements, namely, UAVs, for improving emergency response. To capture a more accurate picture of the problem, a congested model should be considered which adds to the complexity of the problem. Although with the advancement of technology, drones’ characteristics such as endurance and speed are improving rapidly, making aerial humanitarian relief systems more practical and more economic.

In this paper, an edge-based stochastic immobile FLP is studied wherein the best location and number of the facilities are determined with the objective of minimizing the total costs, travel distances and lost demands. The population was assumed to be uniformly distributed along network edges. Since the problem is proved to be NP-hard, a non-dominated sorting GA was developed to solve the problem. The proposed model was applied to the case study of Tehran and the Pareto solutions resulted from the algorithm was presented which could be beneficial to decision makers in this area. The model could be altered to account for multi-service FLP. The integration of an objective function for minimization of waiting times is recommended for future studies.