# Competitive Relief Supply Chain under the Uncertain Conditions

Authors: Amirabbass Boroumand; Majid Nojavan; Davood Mohammaditabar; Reza Ghaemi
Journal: Discrete Dynamics in Nature and Society
 (2023)

---

## Abstract

Relief operations and planning to implement them are essential due to the unpredictability of natural disaster occurrences and their concomitant damages. The Crisis and Disaster Management Organization’s responsibility is to execute plans, coordinate, control relief, and rescue operations to reduce the impacts and consequences of these disasters. Thus, this organization should attend to the struck regions’ needs with their maximum power as soon as possible. This study has initially reviewed the latest available literature regarding the humanitarian supply chain to determine the research gap. Accordingly, in this research, a mathematical model with a leader-follower approach for relief delivery in the crisis response phase was designed, which took into account the policies of the country of Iran for the distribution of relief items. The amount of inventory of suppliers was also considered as uncertain. Then, to validate it, a numerical example was solved by the metaheuristic algorithm MOEA/D. The results indicated that the designed model is valid and the selected algorithm to solve it has an acceptable performance. Finally, some effective parameters were selected in the model and the sensitivity of the model was evaluated based on their changes.


---

## Body

## 1. Introduction

Natural disasters such as earthquakes, floods, hurricanes, and drought embroil different regions of the world every year and often hurt human lives and properties. Moreover, natural disasters are increasing due to factors like population growth and climate change, and the current facilities are insufficient [1]. More comprehensively, natural disasters are associated with potential long-term economic, environmental, financial, human, and social impacts. These unexpected and severe events hurt economic growth [2] and enhance the uncertainty and challenges for organizations [3]. From the management’s perspective, the detrimental impacts of disasters treat societies over the years, and their financial and human costs are inevitable for the victims. However, decision-makers can reduce human recompense and mortalities by their exact, effective, quick, and up-to-date responses. The performance of humanitarian organizations in encountering natural disasters (such as earthquakes, floods, and hurricanes) and human factors (such as fire, environmental pollution, and war) should be assessed in different aspects to optimize the speed of operations, increase the flexibility of the services, and reduce the costs [4]. To this end, disaster management emphasizes planning, prioritization, and decision-making in relief operations [5]. Hence, managing disasters and crises is divided into several steps or phases (commonly included in predisaster and postdisaster phases). Predicting probable properties and lives recompenses and providing preparation programs to reduce the impacts of disasters through improving humanitarian logistics and emergency services are carried out in the predisaster step. In contrast, the postdisaster step comprises measures to respond to the disasters’ recompenses [6]. The mentioned two steps are categorized into four steps for considering more details. The first step is the required measures for preventing a disaster and reducing its disastrous impacts. The second step is preparation measures, including planning manners in society to respond to recompenses quickly. The third step is the response, including utilizing emergency plans and allocated resources to the quick rescue of the victims, providing medical care and treatment services, sending required services and products to injured regions, and helping prevent infrastructural and environmental recompenses. The fourth step, the final recovery, includes measures for returning to the normal situation [7].Based on the discussed topics here, the purpose of this research is to design a mathematical model for planning and organizing relief operation during an earthquake, so as to minimize the costs, harmful environmental effects, and the number of unresponded demands. Based on this, we will first introduce the parameters and decision variables, then we will present the mathematical model, and after that we will solve the model using MOEA/D metaheuristic algorithm, and finally, we will analyze the sensitivity of some parameters.


## 2. Survey on Related Work

Maharjan and Hanaoka [8] introduced a multiobjective location-allocation model for supplying and distributing relief. Their proposed model considers parameters’ inexact and time-varying nature and time-varying envelopment. Results indicated the time, location, and the number of temporary relief centers that should be created. Besides, how the resources should be allocated is revealed in their results. Salehi and Jabarpour [9] assessed a multiobjective model for locating, distributing, and routing multiperiod problems, considering evacuating injured and homeless individuals as well as fuzzy paths in relief services. Notably, the fuzzy approach is used to create uncertain conditions. Jaggi and Singh [10] designed an inventory relief supply chain model for relief product distribution to determine the optimal number of central distribution and relief centers. In this study, a center is considered for relief product distribution. Finally, the model is revised to optimize the total cost of relief operations. It is worth noting that an algorithm is proposed to find the model’s optimal solution. In another study, Zahedi et al. [11] developed two innovative approaches for designing a relief supply chain network using the Internet of Things (IoT) to examine several suspect cases during pandemic (e.g., COVID-19). Their first approach (prioritizing approach) minimizes the maximum response time of ambulances. Moreover, their second approach (allocation approach) minimizes the total crisis response time. Each approach is examined and proved using many experimental problems and one real problem in Iran. Eventually, both of these approaches are proposed. The confirmed cases are reduced by 34.54%, implementing the proposed model based on the IoT in three continuous weeks. Saatchi et al. [12] presented a relief supply chain network in two phases (predisaster and postdisaster). They designed a relief supply chain in two directions (forward and backward directions). Eventually, they solved the model in extensive dimensions using the nondominated sorting genetic algorithm (NSGA-II). Manopiniwes and Irohara [13] studied creating an integrated humanitarian supply chain management model to respond to floods. Their study discusses the interaction among various factors in the relief supply chain during an optimal framework. For this purpose, a model is created to control the total flow of supply distribution, evacuation planning, and relief resources and optimally formulate the routing of the temporary storehouse centers. Finally, a routing model for temporary storehouse problems (created for such disasters) is proposed. Meanwhile, the proposed model is formulated using a multiperiod approach. Mamashli et al. [14] studied the allocation routing problem in the crisis reaction phase. They proposed a scenario-based multiobjective planning model for examining the sustainable allocation of routing problem to cover factors such as sustainability and resilience (rarely considered in previous studies). Moreover, they proposed their model based on the concept of justice, and its purposes are minimizing the travel time, total environmental consequences, and losing total demand. Finally, they proposed a hybrid approach based on a multioptions objective planning method and heuristic solution algorithm to solve the problem in a reasonable time, given the complexity of the model. Vahdani et al. [15] presented a two-objective optimization model to plan a humanitarian regional logistics network. The mentioned model has been considered the extensive range of simultaneous decision-making regarding allocation of emergency facilities location, reclassification, sharing services, and vehicle routing. The two types of open and closed vehicle routing problems are used for ground and aerial routing. Due to the undetermined nature of the disaster, the cost, supply, and demand parameters are considered uncertain in their study. Finally, a hybrid robust optimization model is proposed. The validity of this model is examined in a real study case. Gutierrez et al. [16] assessed the routing problem using automobile routing problem with cumulative capacity. All related constraints are discussed and presented, given the various features of the problem during the time. Additionally, the analysis of related studies regarding solution algorithms and benchmark samples is provided in this study. Finally, some suggestions are provided for future studies in this regard. Emami et al. [17] debated a new term “secondary disaster” recently introduced regarding the consequences of disasters threatening human lives. All hurt regions may encounter secondary disasters due to natural disasters. However, this fact is generally ignored and can amplify the disasters’ consequences. Therefore, they examined the designing of relief resource supply chain to reduce human and economic consequences of natural disasters and presented optimal relief and rescue operations, given the probable occurrence of primary and secondary disasters. Moreover, a single-objective planning mixed nonlinear integer planning model is presented to meet the demand for relief products, victims rescue, and evacuation of injured people, given the dynamic nature of secondary impacts and requiring continuous updating of the relief management process. The mentioned model minimizes the demand in primary and secondary crisis conditions, transportation time, transportation costs, and unmet demand based on the region’s prioritization. Also, the priority of demand scores is defined based on the unmet needs value and the period of privation from relief products and services. According to this perspective, the demand regions are prioritized, and unmet needs are minimized. A hybrid genetic algorithm (GA) approach and rolling-horizon planning are introduced to solve the model as the problem of the study is related to the classification of NP-hard problems. Eventually, the proposed algorithm is implemented on available datasets based on the case study, which indicates the high quality of the solution method in terms of the quality of solution and computational time. Modiri et al. [18] introduced a multiobjective mixed-integer mathematical planning model for designing relief product distribution networks in disaster relief logistics. The first objective function minimizes the total network costs and is divided into two parts: (1) relief costs (transportation, inventory, and fixed facilities costs) and (2) social costs (privation costs). The second objective function minimizes the produced pollution by the network. According to investigating the relevant literature, they declared that this is the first study that proposes a robust fuzzy optimization approach for the problem of designing the relief products distribution network, considering the environmental (CO2 emissions), social (privation cost), and economic consequences under reliability and uncertainty. The multiobjective model is solved using multioption objective planning. A case study based on real data (flood in Sari Province in 2019) is evaluated to demonstrate the model’s validity. The proposed model allows managers and decision-makers to adopt strategic and tactic decisions with the lowest cost and time. Moreover, they can enhance the structure of distribution networks and inventory and reduce victims’ dissatisfaction. Ozen and Krishnamurthy [19] explained that the distribution of relief products for victims is one of the critical activities in reacting to disasters. Due to the dynamic changes in victims’ requirements and postcrisis conditions, the distribution of relief products is challenging. They have modeled distribution operations. So that products such as tarps and blankets are distributed among the victims in various temporary distribution regions called relief centers by volunteers. They examined the impacts of victims’ mobility on the distribution performance. Thus, they modeled every relief center as a queue and distribution operations as a generalized queue network (Gelenbe network). Meanwhile, the product form solutions are examined for the generalized queue network model; a new product form is proved for generalized queue networks with classified batch signals and transmission under particular conditions. They used this result to develop product form approximations applied to a wide range of settings. Then, they applied a generalized queue network model for a case study from relief distribution data of Nepal earthquakes and determined the impact of victims’ mobility on the network performance. Mahmoodi et al. [20] declared that humanitarian supply chain management differs from commercial supply chains. Thus, the purpose of humanitarian supply chains is to minimize the time of reacting to a disaster. Contrary, the purpose of the commercial supply chain is the maximization of profit. Hence, they developed a relief chain structure, which includes emerging technologies in humanitarian provision in two steps (preparation and reaction steps). They sought to maximize the total demand covered by additive manufacturing and distribution centers and the total real weight allocated to drones. Eventually, the proposed model is solved using three methods, including one exact method and two metaheuristic methods. Based on the implementation results, the performance of the nondominated ranking genetic algorithm is better in finding optimal solutions. After solving the model with Cuckoo optimization algorithm (COA) and comparing the results with GAMS software results, the genetics algorithm is better than other options. Hashemi et al. [21] declared that there are still vital and fundamental challenges to optimizing emergency medical services, despite the extensive research and efforts in this regard over the last four decades. The experienced operational problems throughout the implementation of the proposed model led to novel concerns and the development of other models. Therefore, they examined the optimal location of emergency medical centers to present quicker and more efficient care. They sought to propose a mathematical model for the location of medical emergency centers aiming to increase the quality and quantity of demand coverage. Then, the behavior of the model is analyzed by defining numerical examples. The model is solved with GAMS software in small dimensions. The model is calculated in large dimensions using a metaheuristic algorithm by genetics algorithm, based on its NP-hard nature. Eventually, the results of the figures are compared with each other. According to the results, GAMS software’s ability to solve the problem is lost with enhancing the dimensions of the problem; the time of solving the problem is reduced using genetics algorithms compared to GAMS software. Finally, contour lines are used for data analysis in a numerical example. The potential points for emergency medical services followed these lines and performed as the demand points. The accuracy of the model is proved with different parameters. Consequently, their proposed model can respond to medical emergency service demands and determine the optimal location of medical emergency care facilities. More details about survey on related work are shown in Table 1.Table 1
Survey on related work.
Author(s)YearStudied scopeStudied phaseType of goalProblem featuresSolving methodLocating,routing, andtransportationCompetitivenessUncertaintyPreparationPreventionResponse(reaction)ReconstructionEconomicSocialEnvironmentalTimeminimizationMultiproductsMultipletransportationsystemMultiperiodSingle-objectiveMultiobjectiveExactHeuristicMaharjan&Hanaka2020∗∗∗∗∗∗∗Salehi &Jabbarpour2020∗∗∗∗∗∗∗Jaggi& Singh2020∗∗∗∗∗Zahedi et al.2021∗∗∗∗∗Madani-Saatchi et al.2021∗∗∗∗∗∗∗Manopiniwes&Irohara2021∗∗∗∗∗∗Mamashli et al.2021∗∗∗∗∗∗∗∗∗Vahdani et al.2022∗∗∗∗∗∗∗Corona-Gutierrez et al.2022∗∗∗∗∗∗Emami et al.2022∗∗∗∗∗∗∗Modiri et al.2022∗∗∗∗∗∗∗∗∗Ozen& Krishnamurthy2022∗∗∗∗∗∗Mahmoodi et al.2022∗∗∗∗∗∗∗∗∗∗Hashemi et al.2022∗∗∗∗∗Present study∗∗∗∗∗∗∗∗∗∗∗∗


## 3. Research Gap

Based on the above table and considering the examined scopes of previous studies, the considered scopes of the present study are more extensive and examine the considered problem from different perspectives. Also, in the current research, the humanitarian supply chain has been organized locally and in adherence to the policies of the Crisis Management Organization of Iran; accordingly, all public aid and purchased relief items are handed over to the government, which, in turn, will ultimately decide on the quantity and modality of distributing them among the different involved organizations in the relief process. In addition, based on the level of response to the demand in each phase, nonprofit organizations that have cooperated in the distribution process will be provided with an incentive.


### 3.1. Problem Statement

The occurrence of natural disasters, particularly earthquakes, has uncertain nature. Hence, examining conditions and complete preparation in the reaction phase is very important and needs particular attention. Besides, reducing economic, environmental, and mental consequences is essential, needing deep examination and effective and logical solutions [22]. The proper and effective distribution of products and relief services after a disaster occurrence is critical among the mentioned items [23].Therefore, a relief supply chain with triplet levels is introduced in this phase. According to Figure1, in the first level, the supply centers are created permanently far from the disaster location to reduce costs and detrimental environmental consequences (due to relief measures, which are responsible for supplying and sending different types of relief products to distributors). At the second level, distribution centers (public and private) in hurt locations maximize the response temporarily to injured needs. Eventually, demand points are present as the final receiver of relief products.Figure 1
Levels of the proposed supply chain model.


### 3.2. Assumptions

(1)
The inventory of the suppliers is considered to be uncertain (fuzzy).(2)
The costs of supplying and maintaining the products for suppliers have been overlooked.(3)
The costs of establishing distribution centers have been overlooked given the nature of public and private distributors’ establishment (temporary).(4)
The model is multiproduct.(5)
The model is multiperiod.(6)
Various transportation systems (land, marine, and aerial) have been considered.(7)
The possibility of a secondary disaster has not been considered.(8)
The relief products demanded by the distributors (private and public) can be obtained from several supply centers.(9)
Given the diversity of the transportation system, vehicles are assumed to be available in diverse forms, adequate numbers, and sufficient capacity.(10)
The private distribution centers will receive financial incentives from the suppliers per response level (this is the difference between private distribution centers). Meanwhile, the speed of private distributors (the response time to demand) is more than public distributors regarding response to demand.(11)
The transportation costs and other related costs are similar for private and public distributors, as they provide services together.This study tries to increase the coincidence degree of the model with real conditions considering logical conditions to design the model.Optimizing the time and efficiently providing services are critical in the humanitarian supply chain. Hence, the multiplication of supply chain levels leads to the reduction of efficiency and extends the total time of the process. Thereby, the supply chain set does not perform properly.Notably, relief products have high diversity. Hence, the studied model is designed in terms of multiproducts. According to the information regarding supply centers, these centers are used as gyms, storehouses, etc., in emergency cases. The supply centers can earn incomes by this utilization change; hence, this study ignores the costs of supply centers.The sets and indices are as follows:O the set of relief services represented by o indices; o∈O.P the set of periods represented by p indices; p∈P.I: the set of supply centers represented byi indices; i∈I.J: the set of public and private distribution centers represented byj and j′ indices; j∪j′∈J.K the set of demand centers represented by k indices; k∈K.R: the set of transportation represented byr; r∈R.


### 3.3. Parameters

tijrpois the average delivery time spent on relief service o from supply center i to distribution center j (public) using transport method r over the p period.tij′r′po is the average delivery time on relief service o from supply center i to distribution center (private) j′ using transport method r over the p period.tjkr″po is the average delivery time on relief service o from supply center j to demand center (public) k using transportation method r over the p period.tj′kr″′po is the average delivery time on relief service o from supply center j′ to demand center (private) k using transportation method r over the p period.Gijr is the CO2 emission to deliver relief services from supply center i to distribution center j (public) using transportation method r per hour.Gij′r′ is the CO2 emission to deliver relief services from supply center i to distribution center j′ (private) using transportation method r per hour.Gjkr″ is the CO2 emission to deliver relief services from supply center j to demand center k (public) using transportation method r per hour.Gj′kr″′ is the CO2 emission to deliver relief services from supply center j′ to demand center k (private) using transportation method r per hour.Sijro is the costs of delivering relief service o from supply center i to distribution center j (public) using transportation method r per unit.Sij′r′o is the costs of delivering relief service o from supply center i to distribution center j′ (private) using transportation method r per unit.Sjkr″o is the costs of delivering relief service o from supply center j to demand center k (public) using transportation method r per unit.Sj′kr″′o is the costs of delivering relief service o from supply center j′ to demand center k (private) using transportation method r per unit.Qipo˜ is the inventory of relief service o at supply center i over the uncertain p period.Djpo is the expected value of relief service o at supply center j (public) over the p period.Dj′′po is the expected value of relief service o at supply center j′ (private) over the p period.Dk″po is the expected value of relief service o at demand center k over the p period.wjp is the weighted value of distribution center j (public) over the p period based on damage severity.wj′′p is the weighted value of distribution center jprivate over the p period based on damage severity.wk″p is the weighted value of demand center k over the p period based on damage severity.φjpo is the acceptable level of receiving relief service o in distribution center j (public) over the p period.φj′′po is the acceptable level of receiving relief service o in distribution center j′ (private) over the p period.φk″po is the acceptable level of receiving relief service o in demand center k over the p period 3.δko is the risk coefficient of waste related to disaster per thousand units of relief service o on the environment.π is the risk coefficient of carbon dioxide emissions per kilogram on the environment.Uo is the conversion coefficient of relief service o to waste.Bj′kpo is The incentives received value by distribution center j′ (private) for the percentage demand of converting relief services o over the p period.T is the total available time.


### 3.4. Decision Variables

xijrpo is the real value of relief service o provided by supply center i to distribution center j (public) using transportation method r over the p period.xij′r′po is the real value of relief service o provided by supply center i to distribution center j′ (private) using transportation method r over the p period.xjkr″po is the real value of relief service o provided by supply center i to demand center j (public) using transportation method r over the P period.xj′kr″′po is the real value of relief service o provided by supply center i to demand center j′ (private) using transportation method r over the p period.


## 3.1. Problem Statement

The occurrence of natural disasters, particularly earthquakes, has uncertain nature. Hence, examining conditions and complete preparation in the reaction phase is very important and needs particular attention. Besides, reducing economic, environmental, and mental consequences is essential, needing deep examination and effective and logical solutions [22]. The proper and effective distribution of products and relief services after a disaster occurrence is critical among the mentioned items [23].Therefore, a relief supply chain with triplet levels is introduced in this phase. According to Figure1, in the first level, the supply centers are created permanently far from the disaster location to reduce costs and detrimental environmental consequences (due to relief measures, which are responsible for supplying and sending different types of relief products to distributors). At the second level, distribution centers (public and private) in hurt locations maximize the response temporarily to injured needs. Eventually, demand points are present as the final receiver of relief products.Figure 1
Levels of the proposed supply chain model.


## 3.2. Assumptions

(1)
The inventory of the suppliers is considered to be uncertain (fuzzy).(2)
The costs of supplying and maintaining the products for suppliers have been overlooked.(3)
The costs of establishing distribution centers have been overlooked given the nature of public and private distributors’ establishment (temporary).(4)
The model is multiproduct.(5)
The model is multiperiod.(6)
Various transportation systems (land, marine, and aerial) have been considered.(7)
The possibility of a secondary disaster has not been considered.(8)
The relief products demanded by the distributors (private and public) can be obtained from several supply centers.(9)
Given the diversity of the transportation system, vehicles are assumed to be available in diverse forms, adequate numbers, and sufficient capacity.(10)
The private distribution centers will receive financial incentives from the suppliers per response level (this is the difference between private distribution centers). Meanwhile, the speed of private distributors (the response time to demand) is more than public distributors regarding response to demand.(11)
The transportation costs and other related costs are similar for private and public distributors, as they provide services together.This study tries to increase the coincidence degree of the model with real conditions considering logical conditions to design the model.Optimizing the time and efficiently providing services are critical in the humanitarian supply chain. Hence, the multiplication of supply chain levels leads to the reduction of efficiency and extends the total time of the process. Thereby, the supply chain set does not perform properly.Notably, relief products have high diversity. Hence, the studied model is designed in terms of multiproducts. According to the information regarding supply centers, these centers are used as gyms, storehouses, etc., in emergency cases. The supply centers can earn incomes by this utilization change; hence, this study ignores the costs of supply centers.The sets and indices are as follows:O the set of relief services represented by o indices; o∈O.P the set of periods represented by p indices; p∈P.I: the set of supply centers represented byi indices; i∈I.J: the set of public and private distribution centers represented byj and j′ indices; j∪j′∈J.K the set of demand centers represented by k indices; k∈K.R: the set of transportation represented byr; r∈R.


## 3.3. Parameters

tijrpois the average delivery time spent on relief service o from supply center i to distribution center j (public) using transport method r over the p period.tij′r′po is the average delivery time on relief service o from supply center i to distribution center (private) j′ using transport method r over the p period.tjkr″po is the average delivery time on relief service o from supply center j to demand center (public) k using transportation method r over the p period.tj′kr″′po is the average delivery time on relief service o from supply center j′ to demand center (private) k using transportation method r over the p period.Gijr is the CO2 emission to deliver relief services from supply center i to distribution center j (public) using transportation method r per hour.Gij′r′ is the CO2 emission to deliver relief services from supply center i to distribution center j′ (private) using transportation method r per hour.Gjkr″ is the CO2 emission to deliver relief services from supply center j to demand center k (public) using transportation method r per hour.Gj′kr″′ is the CO2 emission to deliver relief services from supply center j′ to demand center k (private) using transportation method r per hour.Sijro is the costs of delivering relief service o from supply center i to distribution center j (public) using transportation method r per unit.Sij′r′o is the costs of delivering relief service o from supply center i to distribution center j′ (private) using transportation method r per unit.Sjkr″o is the costs of delivering relief service o from supply center j to demand center k (public) using transportation method r per unit.Sj′kr″′o is the costs of delivering relief service o from supply center j′ to demand center k (private) using transportation method r per unit.Qipo˜ is the inventory of relief service o at supply center i over the uncertain p period.Djpo is the expected value of relief service o at supply center j (public) over the p period.Dj′′po is the expected value of relief service o at supply center j′ (private) over the p period.Dk″po is the expected value of relief service o at demand center k over the p period.wjp is the weighted value of distribution center j (public) over the p period based on damage severity.wj′′p is the weighted value of distribution center jprivate over the p period based on damage severity.wk″p is the weighted value of demand center k over the p period based on damage severity.φjpo is the acceptable level of receiving relief service o in distribution center j (public) over the p period.φj′′po is the acceptable level of receiving relief service o in distribution center j′ (private) over the p period.φk″po is the acceptable level of receiving relief service o in demand center k over the p period 3.δko is the risk coefficient of waste related to disaster per thousand units of relief service o on the environment.π is the risk coefficient of carbon dioxide emissions per kilogram on the environment.Uo is the conversion coefficient of relief service o to waste.Bj′kpo is The incentives received value by distribution center j′ (private) for the percentage demand of converting relief services o over the p period.T is the total available time.


## 3.4. Decision Variables

xijrpo is the real value of relief service o provided by supply center i to distribution center j (public) using transportation method r over the p period.xij′r′po is the real value of relief service o provided by supply center i to distribution center j′ (private) using transportation method r over the p period.xjkr″po is the real value of relief service o provided by supply center i to demand center j (public) using transportation method r over the P period.xj′kr″′po is the real value of relief service o provided by supply center i to demand center j′ (private) using transportation method r over the p period.


## 4. Mathematical Model

(1)minxijrpo,xjkr″po,xij′r′po,xj′kr″′po∑o∈O∑p∈P4−∑i∈I∑j∈J∑r∈Rwjp×xijrpoDjpo+∑i∈I∑j′∈J∑r∈Rwj′′p×xij′r′poDj′′po−∑j∈J∑k∈K∑r∈Rwkp×xjkr″poDk″po+∑j′∈J∑k∈K∑r∈Rwk″p×xj′kr″′poDk″po,(2)minxijrpo,xjkr″po,xij′r′po,xj′kr″′pο∑o∈O∑p∈P∑i∈I∑j∈J∑r∈RπGijr×tijrpo×xijrpo+∑o∈O∑p∈P∑i∈I∑j′∈J∑r∈RπGij′r′×tij′r′po×xij′r′po+∑o∈O∑p∈P∑j∈J∑k∈K∑r∈RπGjkr″×tjkr″po×xjkr″po+∑o∈O∑p∈P∑j′∈J∑k∈K∑r∈RπGj′kr″′×tj′kr″′po×xj′kr″′po+∑o∈O∑p∈P∑j∈J∑k∈K∑r∈Rδko×Uo×xjkr″po+∑o∈O∑p∈P∑j′∈J∑k∈K∑r∈Rδko×Uo×xj′kr″′po.(3)minxijrpo,xjkr″po,xij′r′po,xj′kr″′po∑o∈O∑p∈P∑i∈I∑j∈J∑r∈RSijro×xijrpo+∑o∈O∑p∈P∑i∈I∑j′∈J∑r∈RSij′r′o×xij′r′po+∑o∈O∑p∈P∑j∈J∑k∈K∑r∈RSjkr″o×xjkr″po+∑o∈O∑p∈P∑j′∈J∑k∈K∑r∈RSj′kr″′o×xj′kr″′po+∑o∈O∑p∈P∑j′∈J∑k∈K∑r∈RBj′kpo×w″kp×xj′kr″′poDk″po.Subject to:(4)∑j∈J∑r∈Rxijrpo+∑j′∈J∑r∈Rxij′r′po=Qipo˜∀i∈I,p∈P,o∈O,(5)∑i∈I∑r∈Rxijrpo≤Djpo∀j∈J,p∈P,o∈O.(6)∑i∈I∑r∈Rxij′r′po≤Dj′po∀j′∈J,p∈P,o∈O.(7)∑i∈I∑r∈Rxijrpo≥φjpoDjpo∀j∈J,p∈P,o∈O.(8)∑i∈I∑r∈Rxij′r′po≥φj′′poDj′′po∀j′∈J,p∈P,o∈O.(9)∑i∈I∑r∈Rtijrpo+∑i∈I∑r∈Rt′ij′rpo+∑k∈K∑r∈Rt″jkrpo+∑k∈K∑r∈Rt‴j′krpo≤T∀j,j′∈J,p∈P,o∈O,(10)xijrpo,xij′r′poare non−negative integer variables∀i∈I,j,j′∈J,r∈R,p∈P,o∈O.(11)maxxjkr″po,xj′kr″′po∑o∈O∑p∈P∑k∈K∑j∈J∑r∈Rwk″p×xjkr″poDk″po×tjkr″po+∑o∈O∑p∈P∑k∈K∑j′∈J∑r∈Rwk″p×xj′kr″′poDk″po×tj′kr″′po.Subject to:(12)∑k∈K∑r∈Rxjkr″po+∑k∈K∑r∈Rxj′kr″′po=∑i∈I∑r∈Rxijrpo+∑i∈I∑r∈Rxij′r′po∀j,j′∈J,p∈P,o∈O.(13)∑j∈J∑r∈Rxjkr″po+∑j′∈J∑r∈Rxj′kr″′po≤Dk″po∀k∈K,p∈P,o∈O.(14)∑j∈J∑r∈Rxjkr″po+∑j′∈J∑r∈Rxj′kr″′po≥φk″poDk″po∀k∈K,p∈P,o∈O.(15)∑j∈J∑r∈Rtjkr″po≥∑j′∈J∑r∈Rtj′kr″′po∀k∈K,p∈P,o∈O.(16)xjkr″po,xj′kr″′poare non−negative integer variables∀j,j′∈J,k∈K,r∈R,p∈P,o∈O.Equation (1) (as the first objective function) minimizes the weighted rate of unmet demand for all periods. Equation (2) is the second objective function of the model in this section, aiming to reduce the detrimental environmental consequences of CO2 emissions and waste. Equation (3) is introduced to minimize the applied costs to the supply chain. Equation (4) ensures that all suppliers’ inventory is delivered to distributors (public and private) in each period. Equation (5) demonstrates that the demand of public distributors is not completely met in every period. Equation (6) demonstrates that the demand of private distributors is not completely met in every period. Equation (7) establishes fair conditions from public distributors as the minimum required demand. Equation (8) establishes fair conditions from private distributors as the minimum required demand. Equation (9) tries to prevent surpassing each period’s total time from the standard determined time. Equation (10) represents the decision variables of the high-level problem.Equation (11) is considered an objective function of low-level and indicates the maximization of the weighted percentage of met demand value in the entire period. Equation (12) demonstrates the balance in each distribution center (public and private) for each period. Equation (13) determines the unmet demand in each period. Equation (14) determines the fair conditions in each period. Equation (15) demonstrates that the total time of relief products and service distribution by private distribution centers is shorter than public distribution centers in each period. Equation (16) demonstrates the decision variables of the low-level problem.


### 4.1. Problem-Solving Approach

Two-level planning is a specific type of multilevel planning. In these models, a subset of the high-level variables is dependent on the solution of the low-level variables. The model is thus solved as follows.The first step: formulation of the integer, two-level, and fuzzy three-objective programming model.According to equation (4), the problem is under uncertain conditions. Thus, the model is made certain at this step. For this purpose, equation (4) is rewritten as follows:(17)∑j∈J∑r∈Rxijrpo+∑j′∈J∑r∈Rxij′r′po=EQipo˜∀i∈I,p∈P,o∈O,where EQipo˜=QiLpo+2QiHpo+QiUpo/4. QiLpo and QiHpo are the pessimistic and optimistic high and low values of inventory, respectively, and QiUpo is the most probable value [25].The second step is the creation of the low-level dual problem.The dual theory is used at this step to form the low-level dual problem. It must be noted that if variablesxijrpo and xij′r′po are constant, the low-level problem can be considered a common transportation problem. Moreover, equation (16) concerning the decision variables of the low-level problem can be reduced as xjkr″po,xj′kr″′po≥0. For this purpose, bjpo,ckpo,dkpo,ekpo are introduced to represent dual variables from equations (12) to (15). As a result, the dual problem related to low-level is given as follows:(18)minbjpo,ckpo,dkpo,ekpo∑o∈O∑p∈P∑j∈J∑i∈I∑r∈Rxijrpo×bjpo+∑o∈O∑p∈P∑j′∈J∑i∈I∑r∈Rxij′r′po×bj′po+∑o∈O∑p∈P∑k∈KDk″po×ckpo+∑o∈O∑p∈P∑k∈K−φkpo×Dk″po×dkpo+∑o∈O∑p∈P∑j′∈J∑k∈K∑r∈R−tj′kr″′po×ekpo.Subject to:(19)bjpo+ckpo−dkpo−ekpo≥wk″pDk″po×tjkr″po∀j∈J,k∈K,r∈R,p∈P,o∈O.(20)bj′po+ckpo−dkpo−ekpo≥wk″pDk″po×tj′kr″′po∀j′∈J,k∈K,r∈R,p∈P,o∈O.(21)bjpo,bj′pours∀j,j′∈J,p∈P,o∈O.(22)ckpo≥0∀k∈K,p∈P,o∈O.(23)dkpo≥0∀k∈K,p∈P,o∈O.(24)ekpo≥0∀k∈K,p∈P,o∈O.Equation (18) is the objective function of the low-level dual problem. Equations (19) to (24) are the constraints of the low-level dual problem (make the solution space of dual problem).The third step: theM1 model is converted into a single-level model with nonlinear equations (M2).(25)minxijrpo,xjkr″po,xij′r′po,xj′kr″′po,bjpo,bj′po,ckpo,dkpo,ekpo∑o∈O∑p∈P4−∑i∈I∑j∈J∑r∈Rwjp×xijrpoDjpo+∑i∈I∑j′∈J∑r∈Rwj′′p×xij′r′poDj′′po−∑j∈J∑k∈K∑r∈Rwk″p×xjkr″poDk″po+∑j′∈J∑k∈K∑r∈Rwk″p×xj′kr″′poDk″po.(26)minxijrpo,xjkr″po,xij′r′po,xj′kr″′po,bjpo,bj′po,ckpo,dkpo,ekpo∑o∈O∑p∈P∑i∈I∑j∈J∑r∈RπGijr×tijrpo×xijrpo+∑o∈O∑p∈P∑i∈I∑j′∈J∑r∈RπGij′r′×tij′r′po×xij′r′po+∑o∈O∑p∈P∑j∈J∑k∈K∑r∈RπGjkr″×tjkr″po×xjkr″po+∑o∈O∑p∈P∑i∈I∑j′∈J∑r∈RπGj′kr″′×tj′kr″′po×xj′kr″′po+∑o∈O∑p∈P∑j∈J∑k∈K∑r∈Rδko×Uo×xjkr″po+∑o∈O∑p∈P∑j′∈J∑k∈K∑r∈Rδko×Uo×xjkr″′po,(27)minxijrpo,xjkr″po,xij′r′po,xj′kr″′po,bjpo,bj′po,ckpo,dkpo,ekpo∑o∈O∑p∈P∑i∈I∑j∈J∑r∈RSijro×xijrpo+∑o∈O∑p∈P∑i∈I∑j′∈J∑r∈RSij′r′o×xij′r′po+∑o∈O∑p∈P∑j∈J∑k∈K∑r∈RSjkr″o×xjkr″po+∑o∈O∑p∈P∑j′∈J∑k∈K∑r∈RSj′kr″′o×xj′kr″′po+∑o∈O∑p∈P∑j′∈J∑k∈K∑r∈RBj′kpo×wk″p×xj′kr″′poDk″po.Subject to:(28)xjkr″po×bjpo+ckpo−dkpo−ekpo−wk″pDk″po×tjkr″po=0∀j∈J,k∈K,r∈R,p∈P,o∈O.(29)xj′kr″′po×bj′po+ckpo−dkpo−ekpo−wk″pDk″po×tj′kr″′po=0∀j′∈J,k∈K,r∈R,p∈P,o∈O.(30)ckpo×Dkpo−∑j∈J∑r∈Rxjkr″po+∑j′∈J∑r∈Rxj′kr″′po=0∀k∈K,p∈P,o∈O,(31)dkpo×∑j∈J∑r∈Rxjkr″po+∑j′∈J∑r∈Rxj′kr″′po−φk″po×Dk″po=0∀k∈K,p∈P,o∈O.(32)ekpo×∑j∈J∑r∈Rtjkr″po−∑j′∈J∑r∈Rtj′kr″′po=0∀k∈K,p∈P,o∈O.(33)xjkr″po≥0∀j∈J,k∈K,r∈R,p∈P,o∈O.(34)xj′kr″′po≥0∀j′∈J,k∈K,r∈R,p∈P,o∈O.Including equations (5)-(10), (12)-(15), (17), and (19)-(24).Equations (25)-(27) are the objective functions of the model that have transformed into a single-level M2. Further, equations (5)-(10), (12)-(15), (17), (33), and (34) specify the feasible region of the primal problem M0 while equations (19)-(24) specify the feasible region of the dual problem. Equations (28)-(32) create sufficient conditions to obtain the optimal value for the primal-dual problems.As a result, theM1 three-objective two-level integer planning model converts to a nonlinear mathematical optimization model, single-level and three-objective, and M2 with xijrpo,xjkr″po,xij′r′po,xj′kr″′po,bjpo,bj′po,ckpo,dkpo,ekpo variables. It must be mentioned that solving the converted model is still challenging since equations (28)-(32) are nonlinear.The fourth step is strategies of converting nonlinear constraints into linear constraints.This step uses binary auxiliary variables to make equations (28)-(32) linear. The N parameter is defined as a large positive constant and θjkrpo∈0,1 is defined as the auxiliary variable for constraint (28). Considering equations (19) and (33), we will have the following:bjpo+ckpo−dkpo−ekpo−wk″p/Dk″po×tjkr″po≥0 and xjkr″po≥0. Consequently, equations (35) and (36) are added as follows:(35)xjkr″po≤N×1−θjkrpo∀j∈J,k∈K,r∈R,p∈P,o∈O.(36)bjpo+ckpo−dkpo−ekpo−wk″pDk″po×tjkr″po≤N×θjkrpo∀j∈J,k∈K,r∈R,p∈P,o∈O.Similarly, considering equations (20) and (34), we will have the following:bj′po+ckpo−dkpo−ekpo−wk″p/Dk″po×tj′kr″′po≥0 and xj′kr″′po≥0. Equations (37) and (38) are thus added as follows:(37)xj′kr″′po≤N×1−θj′krpo∀j′∈J,k∈K,r∈R,p∈P,o∈O.(38)bj′po+ckpo−dkpo−ekpo−wk″p/Dk″po×tj′kr″′po≤N×θj′krpo∀j′∈J,k∈K,r∈R,p∈P,o∈O.In this section,τkpo∈0,1 is defined to make equation (30) linear based on equations (13) and (22). Thus, we will have the following:Dk″po−∑j∈J∑r∈Rxjkr″po+∑j′∈J∑r∈Rxj′kr″′po≥0 and ckpo≥0. Equations (39) and (40) are thus added as follows:(39)ckpo≤N×1−τkpo∀k∈K,p∈P,o∈O.(40)Dk″po−∑j∈J∑r∈Rxjkr″po+∑j′∈J∑r∈Rxj′kr″′po∑j′∈J∑r∈Rxj′kr″′po≤N×τkpo∀k∈K,p∈P,o∈O.εkpo∈0,1 is used to make equation (31) linear based on equations (14) and (23). Thus, we will have: ∑j∈J∑r∈Rxjkr″po+∑j′∈J∑r∈Rxj′kr″′po−φk″po×Dk″po≥0 and dkpo≥0. Equations (41) and (42) are thus added as follows:(41)dkpo≤N×1−εkpo∀k∈K,p∈P,o∈O.(42)∑j∈J∑r∈Rxjkr″po+∑j′∈J∑r∈Rxj′kr″′po−φk″′po×Dk″po≤N×εkpo∀k∈K,p∈P,o∈O.Eventually,σkpo∈0,1 is used to make equation (32) linear based on equations (15) and (24). Thus, we will have the following:∑j∈J∑k∈K∑r∈Rtjkr″po−∑j′∈J∑k∈K∑r∈Rtj′kr″′po≥0 and ekpo≥0. Equations (43) and (44) are thus added as follows:(43)ekpo≤N×1−σkpo∀k∈K,p∈P,o∈O.(44)∑j∈J∑r∈Rtjkr″po−∑j′∈J∑r∈Rtj′kr″′po≤N×σkpo∀k∈K,p∈P,o∈O.The fifth step is creating a multiperiod integer linear planning model ofM3.(45)minxijrpo,xjkr″po,xij′r′po,xj′kr″′po,bjpo,bj′po,ckpo,dkpo,ekpo,θjkrpo,θj′krpo,τkpo,εkpo,σkpo∑o∈O∑p∈P4−∑i∈I∑j∈J∑r∈Rwjp×xijrpoDjpo+∑i∈I∑j′∈J∑r∈Rwj′′p×xij′r′poDj′′po−∑j∈J∑k∈K∑r∈Rwk″p×xjkr″poDk″po+∑j′∈J∑k∈K∑r∈Rwk″′p×xj′kr″′poDk″pominxijrpo,xjkr″po,xij′r′po,xj′kr″′po,bjpo,bj′po,ckpo,dkpo,ekpo,θjkrpo,θj′krpo,τkpo,εkpo,σkpo∑o∈O∑p∈P∑i∈I∑j∈J∑r∈RπGijr×tijrpo×xijrpo+∑o∈O∑p∈P∑i∈I∑j′∈J∑r∈RπGij′r′×tij′r′po×xij′r′po+∑o∈O∑p∈P∑j∈J∑k∈K∑r∈RπGjkr″×tjkr″po×xjkr″po+∑o∈O∑p∈P∑i∈I∑j′∈J∑r∈RπGj′kr″′×tj′kr″′po×xj′kr″′po+∑o∈O∑p∈P∑j∈J∑k∈K∑r∈Rδko×Uo×xjkr″po+∑o∈O∑p∈P∑j′∈J∑k∈K∑r∈Rδko×Uo×xj′kr″′pominxijrpo,xjkr″po,xij′r′po,xj′kr″′po,bjpo,bj′po,ckpo,dkpo,ekpo,θjkrpo,θj′krpo,τkpo,εkpo,σkpo∑o∈O∑p∈P∑i∈I∑j∈J∑r∈RSijro×xijrpo+∑o∈O∑p∈P∑i∈I∑j′∈J∑r∈RSij′r′o×xij′r′po+∑o∈O∑p∈P∑j∈J∑k∈K∑r∈RSjkr″o×xjkr″po+∑o∈O∑p∈P∑j′∈J∑k∈K∑r∈RSj′kr″′o×xj′kr″′po+∑o∈O∑p∈P∑j′∈J∑k∈K∑r∈RBj′kpo×wk″p×xj′kr″′poDk″po.Subject to:(46)∑j∈J∑r∈Rxijrpo+∑j′∈J∑r∈Rxij′r′po=QiLpo+2QiHpo+QiUpo4∀i∈I,p∈P,o∈O.(47)∑i∈I∑r∈Rxijrpo≤Djpo∀j∈J,p∈P,o∈O.(48)∑i∈I∑r∈Rxij′r′po≤Dj′′po∀j′∈J,p∈P,o∈O.(49)∑i∈I∑r∈Rxijrpo≥φjpoDjpo∀j∈J,p∈P,o∈O.(50)∑i∈I∑r∈Rxij′r′po≥φj′′poDj′′po∀j′∈J,p∈P,o∈O.(51)∑i∈I∑r∈Rtijrpo+∑i∈I∑r∈Rtij′r′po+∑k∈K∑r∈Rtjkr″po+∑k∈K∑r∈Rtj′kr″′po≤T∀j,j′∈J,p∈P,o∈O,(52)xijrpo,xij′r′poare non−negative integer variables∀i∈I,j,j′∈J,r∈R,p∈P,o∈O,(53)∑k∈K∑r∈Rxjkr″po+∑k∈K∑r∈Rxj′kr″po=∑i∈I∑r∈Rxijrpo+∑i∈I∑r∈Rxij′r″po∀j,j′∈J,p∈P,o∈O,(54)∑j∈J∑r∈Rxjkr″po+∑j′∈J∑r∈Rxj′kr″′po≤Dk″po∀k∈K,p∈P,o∈O,(55)∑j∈J∑r∈Rxjkr″po+∑j′∈J∑r∈Rxj′kr″po≥φk″poDk″po∀k∈K,p∈P,o∈O,(56)∑j∈J∑r∈Rtjkr″po≥∑j′∈J∑r∈Rtj′kr″′po∀k∈K,p∈P,o∈O,(57)xjkr″po,xj′kr″′poare non−negative integer variables∀j,j′∈J,k∈K,r∈R,p∈P,o∈O,(58)bjpo+ckpo−dkpo−ekpo≥wk″pDk″po×tjkr″po∀j∈J,k∈K,r∈R,p∈P,o∈O,(59)bj′po+ckpo−dkpo−ekpo≥wk″pDk″po×tj′kr″′po∀j′∈J,k∈K,r∈R,p∈P,o∈O,(60)bjpo,bj′pours∀j,j′∈J,p∈P,o∈O,(61)ckpo≥0∀k∈K,p∈P,o∈O,(62)dkpo≥0∀k∈K,p∈P,o∈O,(63)ekpo≥0∀k∈K,p∈P,o∈O,(64)xjkr″po≤N×1−θjkrpo∀j∈J,k∈K,r∈R,p∈P,o∈O,(65)bjpo+ckpo−dkpo−ekpo−wk″pDk″po×tjkr″po≤N×θjkrpo∀j∈J,k∈K,r∈R,p∈P,o∈O,(66)xj′kr″′po≤N×1−θj′krpo∀j′∈J,k∈K,r∈R,p∈P,o∈O,(67)bj′po+ckpo−dkpo−ekpo−wk″pDk″po×tj′kr″′po≤N×θj′krpo∀j′∈J,k∈K,r∈R,p∈P,o∈O,(68)ckpo≤N×1−τkpo∀k∈K,p∈P,o∈O,(69)Dk″po−∑j∈J∑r∈Rxjkr″po+∑j′∈J∑r∈Rxj′kr″′po≤N×τkpo∀k∈K,p∈P,o∈O,(70)dkpo≤N×1−εkpo∀k∈K,p∈P,o∈O,(71)∑j∈J∑r∈Rx′′jkrpo+∑j′∈J∑r∈Rx′′′j′krpo−φ′′kpo×D′′kpo≤N×εkpo∀k∈K,p∈P,o∈O,(72)ekpo≤N×1−σkpo∀k∈K,p∈P,o∈O,(73)∑j∈J∑r∈Rtjkr″po−∑j′∈J∑r∈Rtj′kr″′po≤N×σkpo∀k∈K,p∈P,o∈O,(74)θjkrpo∈0,1∀j∈J,k∈K,r∈R,p∈P,o∈O,(75)θj′krpo∈0,1∀j′∈J,k∈K,r∈R,p∈P,o∈O,(76)τkpo∈0,1∀k∈K,p∈P,o∈O,(77)εkpo∈0,1∀k∈K,p∈P,o∈O,(78)σkpo∈0,1∀k∈K,p∈P,o∈O.The sixth step is solving the model.F1, F2, and F3 demonstrate the three-objective model of the M3, respectively. The branch-and-bound approach is then used to solve the single-objective mixed-integer planning model, using either F1 as the unmet demand rate, F2 as the total environmental damage, or F3 as total costs considered. The final value of each objective function will eventually be demonstrated by F1min, F2min, and F3min.In the seventh step, theM3 model is rewritten as a single-objective model. Equation (79) [24] is used for this purpose. This equation is a weighted linear summation for converting a multiobjective problem to a single-objective problem. Thus, we have the following equation:(79)minF=α1F1−F1minF1min+α2F2−F2minF2min+α3F3−F3minF3min.Now, the model is solved and tested considering the above objective function and equations (46)-(78) so that α1, α2, and α3 are the weights of decision-making factors relevant to the sustainable goals.α1+α2+α3=1 and α1>α2=α3 should also be applied, given that the social sustainability goal was more important than environmental and economic sustainability goals for managers [24].


### 4.2. Parameter Setting

It is noteworthy that the multiobjective evolutionary algorithm based on decomposition (MOEA/D) [26] was implemented using MATLAB to solve the model. Accordingly, the optimal algorithm implementation condition is a subset with a minimum number of unresponded demands, destructive environmental effects, and costs. The considered assumptions are as follows: the relief operation is carried out with four suppliers and three distributors (one of them is a private distributor and two others are public distributors); four demand points for three periods. The model is examined via the mentioned assumptions, and the parameters of the model are valued as follows in Tables 2-31:Table 2
Values of the input parameterwjp.
w110.02w120.07w130.2w210.1w220.012w230.06Table 3
Values of the input parameterwj′′p.
w1′10.11w1′20.27w1′30.08Table 4
Values of the input parameterwk″p.
w1″10.09w1″20.2w1″30.21w2″10.36w2″20.6w2″30.11w3″10.17w3″20.1w3″30.13w4″10.08w4″20.9w4″30.26Table 5
Values of the input parameterDjpo.
D111401679D121624630D131994986D211390624D221602519D2311013412D11220268D12230678D13251591D21218794D22230374D23249749D1134974D1237738D13313013D2134652D2237600D23312667Table 6
Values of the input parameterDj′′po.
D1′11412735D1′21630158D1′311050264D1′1220821D1′2233442D1′3254355D1′135389D1′238015D1′3312782Table 7
Values of the input parameterDk″po.
D1″11109479D1″21164218D1″31273697D2″1188289D2″21132434D2″31220724D3″1198884D3″21148326D3″31247211D4″1156505D4″2184758D4″31141263D1″124414D1″226622D1″3211036D2″124767D2″227151D2″3211919D3″123884D3″225827D3″329711D4″124591D4″226887D4″3211477D1″131191D1″231787D1″332979D2″131059D2″231589D2″332648D3″13926D3″231390D3″332317D4″131235D4″231854D4″333090Table 8
Values of the input parameterGijr.
G1110.295G2110.149G3110.138G1210.252G2210.068G3210.189G1120.314G2120.189G3120.368G1220.347G2220.242G3220.247Table 9
Values of the input parameterGij′r′.
G111′0.173G211′0.012G311′0.141G112′0.255G212′0.396G312′0.187Table 10
Values of the input parameterGjkr″.
G111″0.311G121″0.049G131″0.225G131″0.301G211″0.165G221″0.103G231″0.124G241″0.146G112″0.411G122″0.270G131″0.375G131″0.306G212″0.184G222″0.372G232″0.148G242″0.279Table 11
Values of the input parameterGj′kr‴.
G111‴0.561G121‴0.203G131‴0.227G141‴0.312G112‴0.321G122‴0.299G132‴0.286G142‴0.239Table 12
Values of the input parametertijrpo.
t1111119.32t2111133.12t3111115.18t1112112.6t2112119.8t3112114.4t1113131.2t2113138.4t3113152.8t1111216.5t2111222.5t3111224t111229.9t2112213.5t3112216.2t1113223.1t2113231.5t3113256.7t1111326.4t2111364.8t3111338.4t1112329.7t2112343.2t3112372.9t1113336.3t2113389.1t3113352.8t1211130.36t2211122.08t3211117.9t121219.9t2212121.6t3212111.7t1213133.6t2213126.4t3213157.6t1211240.5t2211219.5t3211227t1212211.7t2212214.4t3212224.3t1213227.3t2213233.6t3213237.8t1211331.2t2211336t3211343.2t1212335.1t2212348.6t3212340.5t1213342.9t2213349.5t3213359.4t1121131.02t2121142.30t3121139.48t1122116.5t2122122.5t3122124t1123152.8t2123167.2t3123176.8t1121250.7t2121242.9t31212109.2t1122220.4t2122214.4t3122226.4t1123239.6t2123256.1t3123259.4t1121343.2t2121357.6t3121372t1122343.2t2122361.2t3122357.6t1123391.2t21233114t3123379.8t1221147.94t2221145.12t3221176.14t1222125.5t2222121t3222140.5t1223181.6t2223172t32231129.6t1221250.6t2221254.6t3221281.8t1222216.8t2222221.6t3222216.8t1223246.2t2223259.1t3223272.6t1221350.4t2221361.2t3221375.6t1222350.4t2222372t3222375.6t1223368.4t2223396.9t32233119.7Table 13
Values of the input parametertij′r′po.
t111′1138.4t211′1132.4t311′1149.2t111′218.4t211′219.3t311′1112.3t111′3140.02t211′3141.4t311′3156.58t111′1227.9t211′1227.9t311′1234.2t111′2246.5t211′2246.5t311′2257t111′3252.08t211′3252.08t311′3263.84t111′1369.3t211′1365.1t311′1375.6t111′2328.8t211′2327.36t311′2315.84t111′3350.4t211′3341.58t311′3334.02t112′1176.8t212′1164.8t312′1198.4t112′1130.6t212′2124.3t312′2135.1t112′3195.88t212′3176.14t312′31109.98t112′1260t212′1257t312′1233t112′2269.3t212′2265.1t312′2275.6t112′3293.06t212′3287.42t312′32101.52t112′13132t212′13105.3t312′13152.1t112′2357.12t212′2345.36t312′2365.52t112′3384.24t212′33113.4t312′33126.36Table 14
Values of the input parametertjkr″po.
t'′111119.75t121″1112.48t131″118.19t141″118.58t111″2140.32t121″2132.76t131″2127.72t141″2125.2t111″3120.25t121″3116.2t131″3125.92t141″3118.63t111″1213.5t121″1217.28t131″1211.34t141″1211.88t111″2226.1t121″2222.62t131″2220.88t141″2217.4t111″3215.75t121″3212.6t131″3220.16t141″3214.49t111″1321t121″1323.52t131″1317.64t141″1321.84t111″2335.4t121″2330.42t131″2328.08t141″2323.4t111″3313.5t121″3310.8t131″3317.28t141″3312.42t211″118.2t221″1110.54t231″1121.3t241″119.45t211″2142.14t221″2131.21t231″2121.54t241″2122.36t211″3115t221″3114.5t231″3122.1t241″3117.6t211″1212.1t221″1211.32t231″129.23t241″1210.6t211″2214.5t221″2212.5t231″2219.54t241″2218.05t211″3213.87t221″3214.41t231″3216.97t241″3213.61t211″1318.5t221″1324.5t231″1315.05t241″1319.02t211″235.48t221″237.54t231″2319.54t241″2312.5t211″3310.4t221″338.5t231″3321.2t241″3310.8t112″1120.91t122″1111.73t132″1111.5t142″116.63t112″2135.1t122″2160.84t132″2154.16t142″2146.64t112″3127.72t122″3121.78t132″3129.7t142″3119.8t112″1272.96t122″1240.32t132″1246.08t142″1232.64t112″2244.1t122″2235.28t132″2229.4t142″2235.28t112″3241.16t122″3232.34t132″3244.1t142″3229.4t112″1322.08t122″1325.92t132″1323.04t142″1318.24t112″2379.2t122″2363.36t132″2320.44t142″2336.45t112″3343.68t122″3334.32t132″3346.8t142″3331.2t212″1118.54t222″1117.12t232″1110.52t242″1115.5t212″2115.47t222″2110.05t232″2115.45t242″2125.25t212″3125.14t222″3122.1t232″3128.4t242″3111.5t212″1212.55t222″1225.25t232″1215.26t242″1229.45t212″2221.45t222″2229.5t232″2221.35t242″2230.48t212″3228.45t222″3227.54t232″3231.21t242″3226.44t212″1324.21t222″1323.87t232″1319.21t242″1316.87t212″2321.15t222″2312.5t232″2315.45t242″2326.44t212″3312.57t222″3310.54t232″3315.84t242″3322.25Table 15
Values of the input parametertj′kr″′po.
t111″′1117t121″′119t131″′119t141″′115t111″′2144t121″′2156t131″′2162t141″′2138t111″′3133t121″′3130.36t131″′3142t141″′3126t111″′126.6t121″′129.6t131″′126.3t141″′127.5t111″′2211t121″′2215.12t131″′2216t141″′2210.26t111″′3233t121″′3233t131″′3239t141″′3226t111″′1335t121″′1320.7t131″′1321.6t141″′1312.6t111″′2314t121″′2318t131″′2320t141″′2312t111″′3323t121″′3325t131″′3332t141″′3320t112″′1153t122″′1128t132″′1134t142″′1120t112″′2130t122″′2138t132″′2142t142″′2126t112″′3150.16t122″′3163t132″′3170t142″′3173t112″′1218t122″′1226t132″′1217t142″′1221t112″′2227.72t122″′2235t132″′2239t142″′2223t112″′3260t122″′3283t132″′3286t142″′3257t112″′1384t122″′1349t132″′1351t142″′1330t112″′2358t122″′2373t132″′2381t142″′2350t112″′3347t122″′3355t132″′3357t142″′3337Table 16
Values of the input parameterδko.
δ110.306δ210.583δ310.253δ410.324δ120.0871δ220.0554δ320.152δ420.234δ130.345δ230.144δ330.1δ430.412Table 17
Values of the input parameterUo.
U10.42U20.25U30.36Table 18
Values of the input parameterSijro.
S1111200000S2111180000S3111200000S1112140000S2112150000S3112140000S11131000000S21131100000S31131200000S1211200000S2211180000S3211200000S1212140000S2212150000S3212120000S12131000000S22132200000S32132000000S112180000S212170000S312180000S112240000S212230000S312240000S1123700000S2123800000S3123700000S122180000S222170000S322180000S122240000S222230000S322240000S1223700000S2223800000S3223700000Table 19
Values of the input parameterSij′r′o.
S111′1250000S211′1200000S311′1250000S111′2150000S211′2140000S311′2150000S111′32000000S211′315000000S311′32000000S112′170000S212′160000S312′170000S112′240000S212′230000S312′240000S112′3800000S212′3700000S312′3800000Table 20
Values of the input parameterSjkr″o.
S111″130000S112″120000S211″140000S212″140000S111″225000S112″217000S211″220000S212″215000S111″3120000S112″380000S211″3140000S212″380000S121″130000S122″120000S221″140000S222″140000S121″225000S122″217000S221″220000S222″215000S121″3120000S122″3100000S221″3140000S222″380000S131″130000S132″120000S231″140000S232″140000S131″225000S132″217000S231″220000S232″215000S131″3120000S132″380000S231″3140000S232″380000S141″130000S142″120000S241″140000S242″140000S141″225000S142″217000S241″220000S242″215000S141″3120000S142″380000S241″3140000S242″380000Table 21
Values of the input parameterSj′kr″′o.
S111″′145000S112″′125000S111″′230000S112″′220000S111″′3125000S112″′385000S121″′145000S122″′125000S121″′230000S122″′220000S121″′3125000S122″′385000S131″′135000S132″′130000S131″′225000S132″′225000S131″′3145000S132″′390000S141″′145000S142″′130000S141″′225000S142″′225000S141″′3145000S142″′390000Table 22
Values of the input parameterBj′kpo.
B11110.4106B11210.7835B11310.2370B11120.2370B11220.5663B11320.5204B11130.3944B11230.3810B11330.4445B12110.3457B12210.6309B12310.6495B12120.6495B12220.2831B12320.3609B12130.7154B12230.5128B12330.3395B13110.3205B13210.3683B13310.0415B13120.0415B13220.8293B13320.8749B13130.2746B13230.3784B13330.7580B14110.6274B14210.1692B14310.1997B14120.1997B14220.3B14320.2655B14130.7265B14230.9051B14330.4120Table 23
Values of the input parameterQiLpo,QiHpo,QiUpo.
Q1L11412848Q1L1216451Q1L135538Q1L21682154Q1L2229455Q1L236837Q1L311083510Q1L3252398Q1L3312412Q1H11522060Q1H1226103Q1H136525Q1H21783091Q1H2239154Q1H239788Q1H311305152Q1H3265257Q1H3316314Q1U11612542Q1U1235761Q1U138025Q1U21894840Q1U2248468Q1U2312235Q1U311487821Q1U3271245Q1U3321371Table 24
Values of the input parameterφjpo.
φ1110.22φ1210.3φ1310.45φ1120.17φ1220.09φ1320.07φ1130.19φ1230.14φ1330.18φ2110.31φ2210.19φ2310.31φ2120.15φ2220.07φ2320.12φ2130.08φ2230.16φ2330.09Table 25
Values of the input parameterφj′′po.
φ1′110.29φ1′210.19φ1′310.18φ1′120.29φ1′220.19φ1′320.18φ1′130.21φ1′230.39φ1′330.14Table 26
Values of the input parameterφk″po.
φ1″110.61φ2″110.51φ3″110.66φ4″110.71φ1″120.09φ2″120.08φ3″120.07φ4″120.15φ1″130.12φ2″130.28φ3″130.19φ4″130.37φ1″210.51φ2″210.47φ3″210.62φ4″210.37φ1″220.28φ2″220.19φ3″220.09φ4″220.17φ1″230.18φ2″230.29φ3″230.37φ4″230.49φ1″310.38φ2″310.61φ3″310.52φ4″310.41φ1″320.08φ2″320.19φ3″320.28φ4″320.36φ1″330.67φ2″330.41φ3″330.59φ4″330.44Table 27
Values of the input parameterbjpo.
b1110.22b1210.11b1310.91b1120.35b1220.34b1320.87b1130.24b1230.08b1330.74b2110.53b2210.37b2310.33b2120.78b2220.87b2320.21b2130.87b2230.14b2330.71Table 28
Values of the input parameterbj′po.
b1110.28b1210.16b1310.4b1120.42b1220.37b1320.32b1130.072b1230.22b1330.61Table 29
Values of the input parameterckpo.
c1110.09c2110.13c3110.45c4110.17c1120.8c2120.02c3120.06c4120.3c1130.81c2130.22c3130.48c4130.24c1210.21c2210.36c3210.21c4210.19c1220.11c2220.75c3220.31c4220.22c1230.08c2230.32c3230.41c4230.10c1310.13c2310.14c3310.51c4310.42c1320.11c2320.62c3320.34c4320.14c1330.52c2330.38c3330.21c4330.12Table 30
Values of the input parameterdkpo.
d1110.56d2110.61d3110.67d4110.23d1120.38d2120.01d3120.58d4120.07d1130.03d2130.64d3130.13d4130.29d1210.28d2210.56d3210.57d4210.92d1220.94d2220.64d3220.60d4220.60d1230.56d2230.12d3230.54d4230.16d1310.27d2310.05d3310.61d4310.72d1320.2d2320.05d3320.44d4320.35d1330.28d2330.32d3330.23d4330.56Table 31
Values of the input parameterekpo.
e1110.75e2110.55e3110.25e4110.04e1120.05e2120.45e3120.36e4120.28e1130.41e2130.5e3130.5e4130.61e1210.5e2210.32e3210.8e4210.67e1220.97e2220.83e3220.5e4220.38e1230.19e2230.11e3230.36e4230.36e1310.24e2310.75e3310.46e4310.27e1320.64e2320.67e3320.05e4320.28e1330.2e2330.57e3330.58e4330.24Further,T = 14400h, π=0.4and r=2.


## 4.1. Problem-Solving Approach

Two-level planning is a specific type of multilevel planning. In these models, a subset of the high-level variables is dependent on the solution of the low-level variables. The model is thus solved as follows.The first step: formulation of the integer, two-level, and fuzzy three-objective programming model.According to equation (4), the problem is under uncertain conditions. Thus, the model is made certain at this step. For this purpose, equation (4) is rewritten as follows:(17)∑j∈J∑r∈Rxijrpo+∑j′∈J∑r∈Rxij′r′po=EQipo˜∀i∈I,p∈P,o∈O,where EQipo˜=QiLpo+2QiHpo+QiUpo/4. QiLpo and QiHpo are the pessimistic and optimistic high and low values of inventory, respectively, and QiUpo is the most probable value [25].The second step is the creation of the low-level dual problem.The dual theory is used at this step to form the low-level dual problem. It must be noted that if variablesxijrpo and xij′r′po are constant, the low-level problem can be considered a common transportation problem. Moreover, equation (16) concerning the decision variables of the low-level problem can be reduced as xjkr″po,xj′kr″′po≥0. For this purpose, bjpo,ckpo,dkpo,ekpo are introduced to represent dual variables from equations (12) to (15). As a result, the dual problem related to low-level is given as follows:(18)minbjpo,ckpo,dkpo,ekpo∑o∈O∑p∈P∑j∈J∑i∈I∑r∈Rxijrpo×bjpo+∑o∈O∑p∈P∑j′∈J∑i∈I∑r∈Rxij′r′po×bj′po+∑o∈O∑p∈P∑k∈KDk″po×ckpo+∑o∈O∑p∈P∑k∈K−φkpo×Dk″po×dkpo+∑o∈O∑p∈P∑j′∈J∑k∈K∑r∈R−tj′kr″′po×ekpo.Subject to:(19)bjpo+ckpo−dkpo−ekpo≥wk″pDk″po×tjkr″po∀j∈J,k∈K,r∈R,p∈P,o∈O.(20)bj′po+ckpo−dkpo−ekpo≥wk″pDk″po×tj′kr″′po∀j′∈J,k∈K,r∈R,p∈P,o∈O.(21)bjpo,bj′pours∀j,j′∈J,p∈P,o∈O.(22)ckpo≥0∀k∈K,p∈P,o∈O.(23)dkpo≥0∀k∈K,p∈P,o∈O.(24)ekpo≥0∀k∈K,p∈P,o∈O.Equation (18) is the objective function of the low-level dual problem. Equations (19) to (24) are the constraints of the low-level dual problem (make the solution space of dual problem).The third step: theM1 model is converted into a single-level model with nonlinear equations (M2).(25)minxijrpo,xjkr″po,xij′r′po,xj′kr″′po,bjpo,bj′po,ckpo,dkpo,ekpo∑o∈O∑p∈P4−∑i∈I∑j∈J∑r∈Rwjp×xijrpoDjpo+∑i∈I∑j′∈J∑r∈Rwj′′p×xij′r′poDj′′po−∑j∈J∑k∈K∑r∈Rwk″p×xjkr″poDk″po+∑j′∈J∑k∈K∑r∈Rwk″p×xj′kr″′poDk″po.(26)minxijrpo,xjkr″po,xij′r′po,xj′kr″′po,bjpo,bj′po,ckpo,dkpo,ekpo∑o∈O∑p∈P∑i∈I∑j∈J∑r∈RπGijr×tijrpo×xijrpo+∑o∈O∑p∈P∑i∈I∑j′∈J∑r∈RπGij′r′×tij′r′po×xij′r′po+∑o∈O∑p∈P∑j∈J∑k∈K∑r∈RπGjkr″×tjkr″po×xjkr″po+∑o∈O∑p∈P∑i∈I∑j′∈J∑r∈RπGj′kr″′×tj′kr″′po×xj′kr″′po+∑o∈O∑p∈P∑j∈J∑k∈K∑r∈Rδko×Uo×xjkr″po+∑o∈O∑p∈P∑j′∈J∑k∈K∑r∈Rδko×Uo×xjkr″′po,(27)minxijrpo,xjkr″po,xij′r′po,xj′kr″′po,bjpo,bj′po,ckpo,dkpo,ekpo∑o∈O∑p∈P∑i∈I∑j∈J∑r∈RSijro×xijrpo+∑o∈O∑p∈P∑i∈I∑j′∈J∑r∈RSij′r′o×xij′r′po+∑o∈O∑p∈P∑j∈J∑k∈K∑r∈RSjkr″o×xjkr″po+∑o∈O∑p∈P∑j′∈J∑k∈K∑r∈RSj′kr″′o×xj′kr″′po+∑o∈O∑p∈P∑j′∈J∑k∈K∑r∈RBj′kpo×wk″p×xj′kr″′poDk″po.Subject to:(28)xjkr″po×bjpo+ckpo−dkpo−ekpo−wk″pDk″po×tjkr″po=0∀j∈J,k∈K,r∈R,p∈P,o∈O.(29)xj′kr″′po×bj′po+ckpo−dkpo−ekpo−wk″pDk″po×tj′kr″′po=0∀j′∈J,k∈K,r∈R,p∈P,o∈O.(30)ckpo×Dkpo−∑j∈J∑r∈Rxjkr″po+∑j′∈J∑r∈Rxj′kr″′po=0∀k∈K,p∈P,o∈O,(31)dkpo×∑j∈J∑r∈Rxjkr″po+∑j′∈J∑r∈Rxj′kr″′po−φk″po×Dk″po=0∀k∈K,p∈P,o∈O.(32)ekpo×∑j∈J∑r∈Rtjkr″po−∑j′∈J∑r∈Rtj′kr″′po=0∀k∈K,p∈P,o∈O.(33)xjkr″po≥0∀j∈J,k∈K,r∈R,p∈P,o∈O.(34)xj′kr″′po≥0∀j′∈J,k∈K,r∈R,p∈P,o∈O.Including equations (5)-(10), (12)-(15), (17), and (19)-(24).Equations (25)-(27) are the objective functions of the model that have transformed into a single-level M2. Further, equations (5)-(10), (12)-(15), (17), (33), and (34) specify the feasible region of the primal problem M0 while equations (19)-(24) specify the feasible region of the dual problem. Equations (28)-(32) create sufficient conditions to obtain the optimal value for the primal-dual problems.As a result, theM1 three-objective two-level integer planning model converts to a nonlinear mathematical optimization model, single-level and three-objective, and M2 with xijrpo,xjkr″po,xij′r′po,xj′kr″′po,bjpo,bj′po,ckpo,dkpo,ekpo variables. It must be mentioned that solving the converted model is still challenging since equations (28)-(32) are nonlinear.The fourth step is strategies of converting nonlinear constraints into linear constraints.This step uses binary auxiliary variables to make equations (28)-(32) linear. The N parameter is defined as a large positive constant and θjkrpo∈0,1 is defined as the auxiliary variable for constraint (28). Considering equations (19) and (33), we will have the following:bjpo+ckpo−dkpo−ekpo−wk″p/Dk″po×tjkr″po≥0 and xjkr″po≥0. Consequently, equations (35) and (36) are added as follows:(35)xjkr″po≤N×1−θjkrpo∀j∈J,k∈K,r∈R,p∈P,o∈O.(36)bjpo+ckpo−dkpo−ekpo−wk″pDk″po×tjkr″po≤N×θjkrpo∀j∈J,k∈K,r∈R,p∈P,o∈O.Similarly, considering equations (20) and (34), we will have the following:bj′po+ckpo−dkpo−ekpo−wk″p/Dk″po×tj′kr″′po≥0 and xj′kr″′po≥0. Equations (37) and (38) are thus added as follows:(37)xj′kr″′po≤N×1−θj′krpo∀j′∈J,k∈K,r∈R,p∈P,o∈O.(38)bj′po+ckpo−dkpo−ekpo−wk″p/Dk″po×tj′kr″′po≤N×θj′krpo∀j′∈J,k∈K,r∈R,p∈P,o∈O.In this section,τkpo∈0,1 is defined to make equation (30) linear based on equations (13) and (22). Thus, we will have the following:Dk″po−∑j∈J∑r∈Rxjkr″po+∑j′∈J∑r∈Rxj′kr″′po≥0 and ckpo≥0. Equations (39) and (40) are thus added as follows:(39)ckpo≤N×1−τkpo∀k∈K,p∈P,o∈O.(40)Dk″po−∑j∈J∑r∈Rxjkr″po+∑j′∈J∑r∈Rxj′kr″′po∑j′∈J∑r∈Rxj′kr″′po≤N×τkpo∀k∈K,p∈P,o∈O.εkpo∈0,1 is used to make equation (31) linear based on equations (14) and (23). Thus, we will have: ∑j∈J∑r∈Rxjkr″po+∑j′∈J∑r∈Rxj′kr″′po−φk″po×Dk″po≥0 and dkpo≥0. Equations (41) and (42) are thus added as follows:(41)dkpo≤N×1−εkpo∀k∈K,p∈P,o∈O.(42)∑j∈J∑r∈Rxjkr″po+∑j′∈J∑r∈Rxj′kr″′po−φk″′po×Dk″po≤N×εkpo∀k∈K,p∈P,o∈O.Eventually,σkpo∈0,1 is used to make equation (32) linear based on equations (15) and (24). Thus, we will have the following:∑j∈J∑k∈K∑r∈Rtjkr″po−∑j′∈J∑k∈K∑r∈Rtj′kr″′po≥0 and ekpo≥0. Equations (43) and (44) are thus added as follows:(43)ekpo≤N×1−σkpo∀k∈K,p∈P,o∈O.(44)∑j∈J∑r∈Rtjkr″po−∑j′∈J∑r∈Rtj′kr″′po≤N×σkpo∀k∈K,p∈P,o∈O.The fifth step is creating a multiperiod integer linear planning model ofM3.(45)minxijrpo,xjkr″po,xij′r′po,xj′kr″′po,bjpo,bj′po,ckpo,dkpo,ekpo,θjkrpo,θj′krpo,τkpo,εkpo,σkpo∑o∈O∑p∈P4−∑i∈I∑j∈J∑r∈Rwjp×xijrpoDjpo+∑i∈I∑j′∈J∑r∈Rwj′′p×xij′r′poDj′′po−∑j∈J∑k∈K∑r∈Rwk″p×xjkr″poDk″po+∑j′∈J∑k∈K∑r∈Rwk″′p×xj′kr″′poDk″pominxijrpo,xjkr″po,xij′r′po,xj′kr″′po,bjpo,bj′po,ckpo,dkpo,ekpo,θjkrpo,θj′krpo,τkpo,εkpo,σkpo∑o∈O∑p∈P∑i∈I∑j∈J∑r∈RπGijr×tijrpo×xijrpo+∑o∈O∑p∈P∑i∈I∑j′∈J∑r∈RπGij′r′×tij′r′po×xij′r′po+∑o∈O∑p∈P∑j∈J∑k∈K∑r∈RπGjkr″×tjkr″po×xjkr″po+∑o∈O∑p∈P∑i∈I∑j′∈J∑r∈RπGj′kr″′×tj′kr″′po×xj′kr″′po+∑o∈O∑p∈P∑j∈J∑k∈K∑r∈Rδko×Uo×xjkr″po+∑o∈O∑p∈P∑j′∈J∑k∈K∑r∈Rδko×Uo×xj′kr″′pominxijrpo,xjkr″po,xij′r′po,xj′kr″′po,bjpo,bj′po,ckpo,dkpo,ekpo,θjkrpo,θj′krpo,τkpo,εkpo,σkpo∑o∈O∑p∈P∑i∈I∑j∈J∑r∈RSijro×xijrpo+∑o∈O∑p∈P∑i∈I∑j′∈J∑r∈RSij′r′o×xij′r′po+∑o∈O∑p∈P∑j∈J∑k∈K∑r∈RSjkr″o×xjkr″po+∑o∈O∑p∈P∑j′∈J∑k∈K∑r∈RSj′kr″′o×xj′kr″′po+∑o∈O∑p∈P∑j′∈J∑k∈K∑r∈RBj′kpo×wk″p×xj′kr″′poDk″po.Subject to:(46)∑j∈J∑r∈Rxijrpo+∑j′∈J∑r∈Rxij′r′po=QiLpo+2QiHpo+QiUpo4∀i∈I,p∈P,o∈O.(47)∑i∈I∑r∈Rxijrpo≤Djpo∀j∈J,p∈P,o∈O.(48)∑i∈I∑r∈Rxij′r′po≤Dj′′po∀j′∈J,p∈P,o∈O.(49)∑i∈I∑r∈Rxijrpo≥φjpoDjpo∀j∈J,p∈P,o∈O.(50)∑i∈I∑r∈Rxij′r′po≥φj′′poDj′′po∀j′∈J,p∈P,o∈O.(51)∑i∈I∑r∈Rtijrpo+∑i∈I∑r∈Rtij′r′po+∑k∈K∑r∈Rtjkr″po+∑k∈K∑r∈Rtj′kr″′po≤T∀j,j′∈J,p∈P,o∈O,(52)xijrpo,xij′r′poare non−negative integer variables∀i∈I,j,j′∈J,r∈R,p∈P,o∈O,(53)∑k∈K∑r∈Rxjkr″po+∑k∈K∑r∈Rxj′kr″po=∑i∈I∑r∈Rxijrpo+∑i∈I∑r∈Rxij′r″po∀j,j′∈J,p∈P,o∈O,(54)∑j∈J∑r∈Rxjkr″po+∑j′∈J∑r∈Rxj′kr″′po≤Dk″po∀k∈K,p∈P,o∈O,(55)∑j∈J∑r∈Rxjkr″po+∑j′∈J∑r∈Rxj′kr″po≥φk″poDk″po∀k∈K,p∈P,o∈O,(56)∑j∈J∑r∈Rtjkr″po≥∑j′∈J∑r∈Rtj′kr″′po∀k∈K,p∈P,o∈O,(57)xjkr″po,xj′kr″′poare non−negative integer variables∀j,j′∈J,k∈K,r∈R,p∈P,o∈O,(58)bjpo+ckpo−dkpo−ekpo≥wk″pDk″po×tjkr″po∀j∈J,k∈K,r∈R,p∈P,o∈O,(59)bj′po+ckpo−dkpo−ekpo≥wk″pDk″po×tj′kr″′po∀j′∈J,k∈K,r∈R,p∈P,o∈O,(60)bjpo,bj′pours∀j,j′∈J,p∈P,o∈O,(61)ckpo≥0∀k∈K,p∈P,o∈O,(62)dkpo≥0∀k∈K,p∈P,o∈O,(63)ekpo≥0∀k∈K,p∈P,o∈O,(64)xjkr″po≤N×1−θjkrpo∀j∈J,k∈K,r∈R,p∈P,o∈O,(65)bjpo+ckpo−dkpo−ekpo−wk″pDk″po×tjkr″po≤N×θjkrpo∀j∈J,k∈K,r∈R,p∈P,o∈O,(66)xj′kr″′po≤N×1−θj′krpo∀j′∈J,k∈K,r∈R,p∈P,o∈O,(67)bj′po+ckpo−dkpo−ekpo−wk″pDk″po×tj′kr″′po≤N×θj′krpo∀j′∈J,k∈K,r∈R,p∈P,o∈O,(68)ckpo≤N×1−τkpo∀k∈K,p∈P,o∈O,(69)Dk″po−∑j∈J∑r∈Rxjkr″po+∑j′∈J∑r∈Rxj′kr″′po≤N×τkpo∀k∈K,p∈P,o∈O,(70)dkpo≤N×1−εkpo∀k∈K,p∈P,o∈O,(71)∑j∈J∑r∈Rx′′jkrpo+∑j′∈J∑r∈Rx′′′j′krpo−φ′′kpo×D′′kpo≤N×εkpo∀k∈K,p∈P,o∈O,(72)ekpo≤N×1−σkpo∀k∈K,p∈P,o∈O,(73)∑j∈J∑r∈Rtjkr″po−∑j′∈J∑r∈Rtj′kr″′po≤N×σkpo∀k∈K,p∈P,o∈O,(74)θjkrpo∈0,1∀j∈J,k∈K,r∈R,p∈P,o∈O,(75)θj′krpo∈0,1∀j′∈J,k∈K,r∈R,p∈P,o∈O,(76)τkpo∈0,1∀k∈K,p∈P,o∈O,(77)εkpo∈0,1∀k∈K,p∈P,o∈O,(78)σkpo∈0,1∀k∈K,p∈P,o∈O.The sixth step is solving the model.F1, F2, and F3 demonstrate the three-objective model of the M3, respectively. The branch-and-bound approach is then used to solve the single-objective mixed-integer planning model, using either F1 as the unmet demand rate, F2 as the total environmental damage, or F3 as total costs considered. The final value of each objective function will eventually be demonstrated by F1min, F2min, and F3min.In the seventh step, theM3 model is rewritten as a single-objective model. Equation (79) [24] is used for this purpose. This equation is a weighted linear summation for converting a multiobjective problem to a single-objective problem. Thus, we have the following equation:(79)minF=α1F1−F1minF1min+α2F2−F2minF2min+α3F3−F3minF3min.Now, the model is solved and tested considering the above objective function and equations (46)-(78) so that α1, α2, and α3 are the weights of decision-making factors relevant to the sustainable goals.α1+α2+α3=1 and α1>α2=α3 should also be applied, given that the social sustainability goal was more important than environmental and economic sustainability goals for managers [24].


## 4.2. Parameter Setting

It is noteworthy that the multiobjective evolutionary algorithm based on decomposition (MOEA/D) [26] was implemented using MATLAB to solve the model. Accordingly, the optimal algorithm implementation condition is a subset with a minimum number of unresponded demands, destructive environmental effects, and costs. The considered assumptions are as follows: the relief operation is carried out with four suppliers and three distributors (one of them is a private distributor and two others are public distributors); four demand points for three periods. The model is examined via the mentioned assumptions, and the parameters of the model are valued as follows in Tables 2-31:Table 2
Values of the input parameterwjp.
w110.02w120.07w130.2w210.1w220.012w230.06Table 3
Values of the input parameterwj′′p.
w1′10.11w1′20.27w1′30.08Table 4
Values of the input parameterwk″p.
w1″10.09w1″20.2w1″30.21w2″10.36w2″20.6w2″30.11w3″10.17w3″20.1w3″30.13w4″10.08w4″20.9w4″30.26Table 5
Values of the input parameterDjpo.
D111401679D121624630D131994986D211390624D221602519D2311013412D11220268D12230678D13251591D21218794D22230374D23249749D1134974D1237738D13313013D2134652D2237600D23312667Table 6
Values of the input parameterDj′′po.
D1′11412735D1′21630158D1′311050264D1′1220821D1′2233442D1′3254355D1′135389D1′238015D1′3312782Table 7
Values of the input parameterDk″po.
D1″11109479D1″21164218D1″31273697D2″1188289D2″21132434D2″31220724D3″1198884D3″21148326D3″31247211D4″1156505D4″2184758D4″31141263D1″124414D1″226622D1″3211036D2″124767D2″227151D2″3211919D3″123884D3″225827D3″329711D4″124591D4″226887D4″3211477D1″131191D1″231787D1″332979D2″131059D2″231589D2″332648D3″13926D3″231390D3″332317D4″131235D4″231854D4″333090Table 8
Values of the input parameterGijr.
G1110.295G2110.149G3110.138G1210.252G2210.068G3210.189G1120.314G2120.189G3120.368G1220.347G2220.242G3220.247Table 9
Values of the input parameterGij′r′.
G111′0.173G211′0.012G311′0.141G112′0.255G212′0.396G312′0.187Table 10
Values of the input parameterGjkr″.
G111″0.311G121″0.049G131″0.225G131″0.301G211″0.165G221″0.103G231″0.124G241″0.146G112″0.411G122″0.270G131″0.375G131″0.306G212″0.184G222″0.372G232″0.148G242″0.279Table 11
Values of the input parameterGj′kr‴.
G111‴0.561G121‴0.203G131‴0.227G141‴0.312G112‴0.321G122‴0.299G132‴0.286G142‴0.239Table 12
Values of the input parametertijrpo.
t1111119.32t2111133.12t3111115.18t1112112.6t2112119.8t3112114.4t1113131.2t2113138.4t3113152.8t1111216.5t2111222.5t3111224t111229.9t2112213.5t3112216.2t1113223.1t2113231.5t3113256.7t1111326.4t2111364.8t3111338.4t1112329.7t2112343.2t3112372.9t1113336.3t2113389.1t3113352.8t1211130.36t2211122.08t3211117.9t121219.9t2212121.6t3212111.7t1213133.6t2213126.4t3213157.6t1211240.5t2211219.5t3211227t1212211.7t2212214.4t3212224.3t1213227.3t2213233.6t3213237.8t1211331.2t2211336t3211343.2t1212335.1t2212348.6t3212340.5t1213342.9t2213349.5t3213359.4t1121131.02t2121142.30t3121139.48t1122116.5t2122122.5t3122124t1123152.8t2123167.2t3123176.8t1121250.7t2121242.9t31212109.2t1122220.4t2122214.4t3122226.4t1123239.6t2123256.1t3123259.4t1121343.2t2121357.6t3121372t1122343.2t2122361.2t3122357.6t1123391.2t21233114t3123379.8t1221147.94t2221145.12t3221176.14t1222125.5t2222121t3222140.5t1223181.6t2223172t32231129.6t1221250.6t2221254.6t3221281.8t1222216.8t2222221.6t3222216.8t1223246.2t2223259.1t3223272.6t1221350.4t2221361.2t3221375.6t1222350.4t2222372t3222375.6t1223368.4t2223396.9t32233119.7Table 13
Values of the input parametertij′r′po.
t111′1138.4t211′1132.4t311′1149.2t111′218.4t211′219.3t311′1112.3t111′3140.02t211′3141.4t311′3156.58t111′1227.9t211′1227.9t311′1234.2t111′2246.5t211′2246.5t311′2257t111′3252.08t211′3252.08t311′3263.84t111′1369.3t211′1365.1t311′1375.6t111′2328.8t211′2327.36t311′2315.84t111′3350.4t211′3341.58t311′3334.02t112′1176.8t212′1164.8t312′1198.4t112′1130.6t212′2124.3t312′2135.1t112′3195.88t212′3176.14t312′31109.98t112′1260t212′1257t312′1233t112′2269.3t212′2265.1t312′2275.6t112′3293.06t212′3287.42t312′32101.52t112′13132t212′13105.3t312′13152.1t112′2357.12t212′2345.36t312′2365.52t112′3384.24t212′33113.4t312′33126.36Table 14
Values of the input parametertjkr″po.
t'′111119.75t121″1112.48t131″118.19t141″118.58t111″2140.32t121″2132.76t131″2127.72t141″2125.2t111″3120.25t121″3116.2t131″3125.92t141″3118.63t111″1213.5t121″1217.28t131″1211.34t141″1211.88t111″2226.1t121″2222.62t131″2220.88t141″2217.4t111″3215.75t121″3212.6t131″3220.16t141″3214.49t111″1321t121″1323.52t131″1317.64t141″1321.84t111″2335.4t121″2330.42t131″2328.08t141″2323.4t111″3313.5t121″3310.8t131″3317.28t141″3312.42t211″118.2t221″1110.54t231″1121.3t241″119.45t211″2142.14t221″2131.21t231″2121.54t241″2122.36t211″3115t221″3114.5t231″3122.1t241″3117.6t211″1212.1t221″1211.32t231″129.23t241″1210.6t211″2214.5t221″2212.5t231″2219.54t241″2218.05t211″3213.87t221″3214.41t231″3216.97t241″3213.61t211″1318.5t221″1324.5t231″1315.05t241″1319.02t211″235.48t221″237.54t231″2319.54t241″2312.5t211″3310.4t221″338.5t231″3321.2t241″3310.8t112″1120.91t122″1111.73t132″1111.5t142″116.63t112″2135.1t122″2160.84t132″2154.16t142″2146.64t112″3127.72t122″3121.78t132″3129.7t142″3119.8t112″1272.96t122″1240.32t132″1246.08t142″1232.64t112″2244.1t122″2235.28t132″2229.4t142″2235.28t112″3241.16t122″3232.34t132″3244.1t142″3229.4t112″1322.08t122″1325.92t132″1323.04t142″1318.24t112″2379.2t122″2363.36t132″2320.44t142″2336.45t112″3343.68t122″3334.32t132″3346.8t142″3331.2t212″1118.54t222″1117.12t232″1110.52t242″1115.5t212″2115.47t222″2110.05t232″2115.45t242″2125.25t212″3125.14t222″3122.1t232″3128.4t242″3111.5t212″1212.55t222″1225.25t232″1215.26t242″1229.45t212″2221.45t222″2229.5t232″2221.35t242″2230.48t212″3228.45t222″3227.54t232″3231.21t242″3226.44t212″1324.21t222″1323.87t232″1319.21t242″1316.87t212″2321.15t222″2312.5t232″2315.45t242″2326.44t212″3312.57t222″3310.54t232″3315.84t242″3322.25Table 15
Values of the input parametertj′kr″′po.
t111″′1117t121″′119t131″′119t141″′115t111″′2144t121″′2156t131″′2162t141″′2138t111″′3133t121″′3130.36t131″′3142t141″′3126t111″′126.6t121″′129.6t131″′126.3t141″′127.5t111″′2211t121″′2215.12t131″′2216t141″′2210.26t111″′3233t121″′3233t131″′3239t141″′3226t111″′1335t121″′1320.7t131″′1321.6t141″′1312.6t111″′2314t121″′2318t131″′2320t141″′2312t111″′3323t121″′3325t131″′3332t141″′3320t112″′1153t122″′1128t132″′1134t142″′1120t112″′2130t122″′2138t132″′2142t142″′2126t112″′3150.16t122″′3163t132″′3170t142″′3173t112″′1218t122″′1226t132″′1217t142″′1221t112″′2227.72t122″′2235t132″′2239t142″′2223t112″′3260t122″′3283t132″′3286t142″′3257t112″′1384t122″′1349t132″′1351t142″′1330t112″′2358t122″′2373t132″′2381t142″′2350t112″′3347t122″′3355t132″′3357t142″′3337Table 16
Values of the input parameterδko.
δ110.306δ210.583δ310.253δ410.324δ120.0871δ220.0554δ320.152δ420.234δ130.345δ230.144δ330.1δ430.412Table 17
Values of the input parameterUo.
U10.42U20.25U30.36Table 18
Values of the input parameterSijro.
S1111200000S2111180000S3111200000S1112140000S2112150000S3112140000S11131000000S21131100000S31131200000S1211200000S2211180000S3211200000S1212140000S2212150000S3212120000S12131000000S22132200000S32132000000S112180000S212170000S312180000S112240000S212230000S312240000S1123700000S2123800000S3123700000S122180000S222170000S322180000S122240000S222230000S322240000S1223700000S2223800000S3223700000Table 19
Values of the input parameterSij′r′o.
S111′1250000S211′1200000S311′1250000S111′2150000S211′2140000S311′2150000S111′32000000S211′315000000S311′32000000S112′170000S212′160000S312′170000S112′240000S212′230000S312′240000S112′3800000S212′3700000S312′3800000Table 20
Values of the input parameterSjkr″o.
S111″130000S112″120000S211″140000S212″140000S111″225000S112″217000S211″220000S212″215000S111″3120000S112″380000S211″3140000S212″380000S121″130000S122″120000S221″140000S222″140000S121″225000S122″217000S221″220000S222″215000S121″3120000S122″3100000S221″3140000S222″380000S131″130000S132″120000S231″140000S232″140000S131″225000S132″217000S231″220000S232″215000S131″3120000S132″380000S231″3140000S232″380000S141″130000S142″120000S241″140000S242″140000S141″225000S142″217000S241″220000S242″215000S141″3120000S142″380000S241″3140000S242″380000Table 21
Values of the input parameterSj′kr″′o.
S111″′145000S112″′125000S111″′230000S112″′220000S111″′3125000S112″′385000S121″′145000S122″′125000S121″′230000S122″′220000S121″′3125000S122″′385000S131″′135000S132″′130000S131″′225000S132″′225000S131″′3145000S132″′390000S141″′145000S142″′130000S141″′225000S142″′225000S141″′3145000S142″′390000Table 22
Values of the input parameterBj′kpo.
B11110.4106B11210.7835B11310.2370B11120.2370B11220.5663B11320.5204B11130.3944B11230.3810B11330.4445B12110.3457B12210.6309B12310.6495B12120.6495B12220.2831B12320.3609B12130.7154B12230.5128B12330.3395B13110.3205B13210.3683B13310.0415B13120.0415B13220.8293B13320.8749B13130.2746B13230.3784B13330.7580B14110.6274B14210.1692B14310.1997B14120.1997B14220.3B14320.2655B14130.7265B14230.9051B14330.4120Table 23
Values of the input parameterQiLpo,QiHpo,QiUpo.
Q1L11412848Q1L1216451Q1L135538Q1L21682154Q1L2229455Q1L236837Q1L311083510Q1L3252398Q1L3312412Q1H11522060Q1H1226103Q1H136525Q1H21783091Q1H2239154Q1H239788Q1H311305152Q1H3265257Q1H3316314Q1U11612542Q1U1235761Q1U138025Q1U21894840Q1U2248468Q1U2312235Q1U311487821Q1U3271245Q1U3321371Table 24
Values of the input parameterφjpo.
φ1110.22φ1210.3φ1310.45φ1120.17φ1220.09φ1320.07φ1130.19φ1230.14φ1330.18φ2110.31φ2210.19φ2310.31φ2120.15φ2220.07φ2320.12φ2130.08φ2230.16φ2330.09Table 25
Values of the input parameterφj′′po.
φ1′110.29φ1′210.19φ1′310.18φ1′120.29φ1′220.19φ1′320.18φ1′130.21φ1′230.39φ1′330.14Table 26
Values of the input parameterφk″po.
φ1″110.61φ2″110.51φ3″110.66φ4″110.71φ1″120.09φ2″120.08φ3″120.07φ4″120.15φ1″130.12φ2″130.28φ3″130.19φ4″130.37φ1″210.51φ2″210.47φ3″210.62φ4″210.37φ1″220.28φ2″220.19φ3″220.09φ4″220.17φ1″230.18φ2″230.29φ3″230.37φ4″230.49φ1″310.38φ2″310.61φ3″310.52φ4″310.41φ1″320.08φ2″320.19φ3″320.28φ4″320.36φ1″330.67φ2″330.41φ3″330.59φ4″330.44Table 27
Values of the input parameterbjpo.
b1110.22b1210.11b1310.91b1120.35b1220.34b1320.87b1130.24b1230.08b1330.74b2110.53b2210.37b2310.33b2120.78b2220.87b2320.21b2130.87b2230.14b2330.71Table 28
Values of the input parameterbj′po.
b1110.28b1210.16b1310.4b1120.42b1220.37b1320.32b1130.072b1230.22b1330.61Table 29
Values of the input parameterckpo.
c1110.09c2110.13c3110.45c4110.17c1120.8c2120.02c3120.06c4120.3c1130.81c2130.22c3130.48c4130.24c1210.21c2210.36c3210.21c4210.19c1220.11c2220.75c3220.31c4220.22c1230.08c2230.32c3230.41c4230.10c1310.13c2310.14c3310.51c4310.42c1320.11c2320.62c3320.34c4320.14c1330.52c2330.38c3330.21c4330.12Table 30
Values of the input parameterdkpo.
d1110.56d2110.61d3110.67d4110.23d1120.38d2120.01d3120.58d4120.07d1130.03d2130.64d3130.13d4130.29d1210.28d2210.56d3210.57d4210.92d1220.94d2220.64d3220.60d4220.60d1230.56d2230.12d3230.54d4230.16d1310.27d2310.05d3310.61d4310.72d1320.2d2320.05d3320.44d4320.35d1330.28d2330.32d3330.23d4330.56Table 31
Values of the input parameterekpo.
e1110.75e2110.55e3110.25e4110.04e1120.05e2120.45e3120.36e4120.28e1130.41e2130.5e3130.5e4130.61e1210.5e2210.32e3210.8e4210.67e1220.97e2220.83e3220.5e4220.38e1230.19e2230.11e3230.36e4230.36e1310.24e2310.75e3310.46e4310.27e1320.64e2320.67e3320.05e4320.28e1330.2e2330.57e3330.58e4330.24Further,T = 14400h, π=0.4and r=2.


## 5. System Hardware Specifications for the Calculations

We used a laptop with 8 MB RAM, a 7-core Intel processor, and 64-bit Windows 10 to program and test the proposed method as well as the methods being compared in MATLAB. In addition, there was no specific limit on storage. A 1TB hard drive with approximately 700 GB of free space was utilized. The 2019a version of MATLAB was also used to conduct the tests on the proposed method and the method in comparison.Table32 illustrates the control parameters used in MOEA/D to solve the proposed chain problem. The population size parameter indicates the number of solutions in the MOEA/D search space, which has been considered 20 through a trial and error test. As mentioned, the function of heuristic optimization methods is derived from generating multiple solutions to a problem in the optimization search space, and then, discovering and brainstorming the solutions to extract the proper solution, that is, the best solution as the response to the problem. The reservoir parameter stores the dominant solutions throughout each iteration. Given that there are three objectives when faced with multiple objectives in heuristic optimization algorithms, we require a reservoir to store the most optimal solution in the form of the Pareto archive. The number of objective parameter depicts the number of objectives in the optimization algorithm: the first objective is the minimization of unresponded demands, the second objective is the minimization of destructive environmental effects, and the third objective is the costs of the supply chain. Also, maximum iterations show the limit of the algorithm implementation to obtain convergence, which is also the termination condition of the algorithm. The algorithm continues until it achieves convergence or maximum iterations.Table 32
Control parameters.
MOEA/D optimizationPopulation size20Reservoir15Number of objectives3Maximum iteration250


### 5.1. The Pseudocode of the Proposed Algorithm

A population of N pointx1,…,xN, where xi is the current solution to the ith subproblem.Objective functionF1,…,FN, where Fi is the F-value of xi, that is, Fi=Fxi.Initialize max iteration, that is, max it.While (t < Max it number of iterations.For every solution in population, that is,. i = 1,…, N.(i)
Initialize position of solution based on step 3-3-1-1-1.(ii)
Compute the Euclidean distances between any two weight vectors and then work out the T and closest weight vectors to each weight vector. For each i = 1,…, N, set B(i) = {i1.….iT} whereλi1,…,λiTare T closest weight vectors toλi.(iii)
Evaluate fitness of solution based on step. 3-3-1-1-2.(iv)
z=z1,…,zmTwhereziis the best value found so far objectivefi(three-objective) and T is the number of the weight vectors in the neighborhood of each weight vector.(v)
An external repository (ER) which is used to store nondominated solutions found during the search.(vi)
Select a solution randomly and perform single point mutation.(vii)
Updating locations/solutions step 3-3-1-1-2 using a leader from ER based on best, the best solution is extracted via previous phase.End for.For each j = 1, …, m.(viii)
Updatez=z1,…,zmTs: for each j = 1, …, m, ifzj<fjy′,then setzj=fjy′.(ix)
Update position so ifj∈Bi.ifmaxλifjy′−Zt∗≤maxλifjxj−Zt∗,then setxj=y′andFj=Fy′.End for.Recalculate from (ER) of all the vectors dominated byFy′.Recalculate fitness of all solutions.Rank the solution/solution and find the current bestx∗.End while.Results: best solution.


## 5.1. The Pseudocode of the Proposed Algorithm

A population of N pointx1,…,xN, where xi is the current solution to the ith subproblem.Objective functionF1,…,FN, where Fi is the F-value of xi, that is, Fi=Fxi.Initialize max iteration, that is, max it.While (t < Max it number of iterations.For every solution in population, that is,. i = 1,…, N.(i)
Initialize position of solution based on step 3-3-1-1-1.(ii)
Compute the Euclidean distances between any two weight vectors and then work out the T and closest weight vectors to each weight vector. For each i = 1,…, N, set B(i) = {i1.….iT} whereλi1,…,λiTare T closest weight vectors toλi.(iii)
Evaluate fitness of solution based on step. 3-3-1-1-2.(iv)
z=z1,…,zmTwhereziis the best value found so far objectivefi(three-objective) and T is the number of the weight vectors in the neighborhood of each weight vector.(v)
An external repository (ER) which is used to store nondominated solutions found during the search.(vi)
Select a solution randomly and perform single point mutation.(vii)
Updating locations/solutions step 3-3-1-1-2 using a leader from ER based on best, the best solution is extracted via previous phase.End for.For each j = 1, …, m.(viii)
Updatez=z1,…,zmTs: for each j = 1, …, m, ifzj<fjy′,then setzj=fjy′.(ix)
Update position so ifj∈Bi.ifmaxλifjy′−Zt∗≤maxλifjxj−Zt∗,then setxj=y′andFj=Fy′.End for.Recalculate from (ER) of all the vectors dominated byFy′.Recalculate fitness of all solutions.Rank the solution/solution and find the current bestx∗.End while.Results: best solution.


## 6. Result

One way to ensure the accuracy of the response is to analyze the performance of the optimization algorithms during each iteration (algorithm implementation). The response of performance accuracy is often normal within [0-1], with 1 denoting the 100% accuracy of the algorithm and 0 indicating inaccuracy in finding the optimal response to the problem. Therefore, the closer the algorithm performance is to 1, the more precise it is, and consequently, the more reliable the responses are. In addition, if an algorithm has an upward trend during each implementation (iteration) in terms of performance accuracy, it indicates that the algorithm is highly capable of solving the problem. Figure2 shows the performance accuracy of the algorithm. As can be observed, during the implementation of the algorithm, 1 to 250 iterations have been determined, and the algorithm has shown an upward trend with regards to accuracy. Also, at roughly the 200th iteration, the algorithm achieves 90% convergence, with no changes in the responses from that point forward. This shows that the algorithm performance is highly accurate; therefore, the findings from the implementation of this algorithm can be trusted to solve the supply chain problem being analyzed.Figure 2
Performance trend of the algorithm.Figure3 illustrates the results from the effects of the algorithm implementation on the first objective, that is, the minimization of the unresponded demands. It can be seen that the range of the changes in each iteration of the algorithm, up to the maximum iteration, is [5.5-8.5]. With an increase in the iteration, the unresponded demand converges at approximately 0.5.Figure 3
Not respond demand (first objective) considering the iteration of the algorithm.Figure4 illustrates the results from the effects of the algorithm implementation on the second objective of the model, that is, the minimization of destructive environmental effects. It can be seen that the range of the changes in each iteration of the algorithm, up to the maximum iteration, is [0.1-0.8]. With an increase in the iteration, the destructive environmental effects converge at approximately 0.1.Figure 4
Destructive environmental effects in the supply chain (the second objective) considering the iteration of the algorithm.Figure5 illustrates the results from the effects of the algorithm implementation on the second objective of the model, that is, the minimization of destructive environmental effects. It can be seen that the range of the changes in each iteration of the algorithm, up to the maximum iteration, is [0.2-1]. During an increase in the iteration, the cost converges at approximately 0.2.Figure 5
Costs (the third objective) considering the iteration of the algorithm.Given that the designed relief supply chain in this study is multiobjective, we need to extract the Pareto archive from the algorithm for the comparison of one objective with other objectives. The Pareto archive shows whether the algorithm sought to achieve optimal responses or not. Considering the objective functions introduced in the humanitarian supply chain model in this study, Figure6 illustrates the Pareto archive of each objective.Figure 6
Pareto archive of the three objectives.Figure6 is a 3D depiction of the Pareto archive of the responses obtained during each iteration of the algorithm for every three objectives. According to the graph, the algorithm has clearly attempted to reduce the unresponded demands, the destructive environmental effects, and the costs to solve the problem. Further, the problem has been solved in a way to minimize each of the objectives during each iteration. The majority of the responses have been compressed near the minimum ranges of each three objectives, and the dispersion of the responses is marginal. The results reveal that, ultimately, the responses for the unresponded demands have converged at roughly 5%, for the destructive environmental effects at 10%, and for the costs at 20%. Figures 7-9 show separate 2D graphs of the Pareto archive of the objective functions.Figure 7
Pareto archive from the first and second objective functions.Figure 8
Pareto archive from the first and third objective functions.Figure 9
Pareto archive from the second and third objective functions.Findings obtained from numerous experiments revealed that the weights of the first objective are greater than the other two objectives with the first objective weighing 0.5α1=0.5, the second 0.2 α1=0.2, and the third 0.3 α1=0.3. Also, the values of the decision variables are shown in Tables 33-36, and the changes in each objective function per 50 iterations are shown in Table 37.Table 33
Values of the decision variablexijrpo.
x1111114938x11112580x11113154x1211116296x12112871x12113103x1121113276x11212382x1121362x122118359x12212313x1221359x1112114417x111221307x11123122x1212122655x12122532x12123183x1122111357x11222501x11223112x122217150x12222410x12223107x1113137063x111321549x11133132x1213133975x121321033x12133199x1123122452x11232491x11233160x1223114136x12232401x12233168x2111110864x21112354x2111386x221119506x22112419x2211374x212116883x21212295x2121335x222117375x22212295x2221350x2112113387x21122871x2112388x2212124715x22122726x22123108x212216309x21222273x2122391x222215888x22222387x2222385x2113124709x21132861x21133118x2213121620x22132918x22133110x2123111641x21232379x21233112x222311273x22232267x22233136x311117469x31112484x3111391x321118827x32112516x3211363x312117867x31212208x3121341x322115408x32212243x3221347x3112111327x31122629x3112374x3212116476x32122774x32123102x312216730x31222387x3122364x322214626x32222319x3222375x3113120076x31132746x3113395x3213116987x32132631x3213381x3123112473x31232312x31233128x322319147x32232379x3223396Table 34
Values of the decision variablexij′r′po.
x111′1142792x111′121573x111′13485x112′1132282x112′121552x112′13357x111′2148602x111′222204x111′23709x112′2150084x112′221451x112′23589x111′3189892x111′324183x111′33680x112′3159420x112′321671x112′33500x211′1133399x211′121284x211′13461x212′1125196x212′121422x212′13311x211′2136748x211′221798x211′23674x212′2143663x212′221249x212′23408x211′3163582x211′323412x211′33561x212′3151802x212′321570x212′33449x311′1128180x311′121284x311′13267x312′1121259x312′121336x312′13247x311′2133191x311′221798x311′23390x312′2134673x312′221330x312′23514x311′3165774x311′323412x311′33459x312′3141137x312′321823x312′33333Table 35
Values of the decision variablexjkr″po.
x111″118487x211″1129258x111″12403x211″12194x111″1371x211″1365x112″1110080x212″1120718x112″12320x212″12240x112″1334x212″1346x111″2116476x211″2115476x111″22726x211″22673x111″23204x211″23189x112″2112618x212″2111618x112″22341x212″22300x112″23160x212″23160x111″3117759x211″3116454x111″32717x211″32700x111″33118x211″33100x112″3111641x212″3111541x112″32312x212″32300x212″33112x212″33112x121″1110864x221″112117x121″12516x221″12315x121″1380x221″1376x122″115654x222″112178x122″12184x222″1284x122″1340x222″1350x121″2113384x221″2113267x121″22629x221″22679x121″23163x221″23178x122″21111777x222″21111756x122″22273x222″22314x122″23128x222″23128x121″3124709x221″3123521x121″32918x211″32935x121″3392x221″33110x122″3112473x222″3111545x122″32334x222″32350x122″33120x222″3380x131″117129x231″1113737x131″12338x231″12328x131″1360x231″1367x132″115654x232″11152x132″12210x232″12310x132″1335x232″1330x131″2111327x231″2110654x131″22580x231″22440x131″23176x231″23131x132″2110095x232″219541x132″22273x232″22200x132″23107x232″23150x131″3115443x231″3113527x131″32660x231″32590x131″3384x231″33105x132″319147x232″318441x132″32245x232″32220x132″3380x232″3380x141″117469x241″112866x141″12354x241″12350x141″1374x241″1382x142″11316x242″118830x142″12147x242″12145x142″1328x242″1330x141″211298x241″219756x141″22484x241″22624x141″23136x241″23181x141″217571x242″217462x142″22223x242″22296x142″23128x242″2385x141″3119304x241″3118452x141″32574x241″32644x141″3373x241″3373x142″318315x242″317641x141″32223x242″32248x142″3388x242″33120Table 36
Values of the decision variablexj′kr″′po.
x111″′1130708x111″′121325x111″′1371x112″′1142792x112″′121681x112″′1340x111″′2136748x111″′221279x111″′23496x112″′2135957x112″′221128x112″′23468x111″′3170159x111″′321486x111″′33544x112″′3147231x112″′321063x112″′33308x121″′1111810x121″′12911x121″′1380x122″′1124005x122″′121034x122″′1335x121″′2126079x121″′221624x121″′23549x122″′2128252x122″′221624x122″′23423x121″′3150427x121″′321783x121″′33391x122″′3128948x122″′321519x122″′33359x131″′1119684x131″′121035x131″′1360x132″′1124005x132″′12991x132″′1334x131″′2133191x131″′221798x131″′23390x132″′2139810x132″′221798x132″′23332x131″′315481x131″′321486x131″′33340x132″′3133519x132″′321013x132″′33372x141″′1116535x141″′12869x141″′1374x142″′1113568x142″′12603x142″′1328x141″′2122522x141″′221102x141″′23337x142″′2124399x142″′221102x142″′23287x141″′3143849x141″′321188x141″′33425x142″′3142661x142″′32861x142″′33243Table 37
Changes in each objective function per 50 iterations.
IterationChanges in thefirst objectivefunction(percentage)Changesin thesecondobjectivefunction(normalizedbetween0 and 1)Changes inthe thirdobjectivefunction(Toman)507.990.908734675000001007.80.702144510000001506.80.453984380000002006.40.351933984000002505.50.26031375000000


## 7. Sensitivity Analysis

To examine the effect of parameterwk″p on the first and third objective functions, we solved the problem with different values (+10%, +5%, ground state, −5%, and −10%) and reported the behavior of the first and third objective functions. According to Figure 10, a 5% increase in the parameter does not have a significant impact on the first objective function. However, a decrease by the same amount would reduce in first objective function by 30%. Further, a 10% increase/decrease in the parameter on the third objective function is substantial.Figure 10
Sensitivity analysis of parameterwk″p and its impacts on the first and third objective functions.Figure11 shows that a 5% decrease/increase in parameter δko does not cause a significant change in the second objective function. However, a 10% decrease/increase in the coefficient of the parameter would decrease and increase the second objective function by 17% and 12%, respectively.Figure 11
Sensitivity analysis of parameterδkp and its impacts on the second objective function.Figure12 shows the sensitivity analysis of parameter Dk″p on the first and third objective functions. The increase/decrease in the perimeter on the third objective function in the considered scope is not substantial, whereas, for instance, a 20% increase in the parameter on the first objective function is.Figure 12
Sensitivity analysis of parameterDk″p and its impacts on the first and third objective functions.


## 8. Conclusions and Outlook

In this research, an attempt was made to provide a design by observing the dimensions of the sustainability in the humanitarian supply chain model in the face of the crisis, in a local way according to the policies of the country of Iran. Considering the research gap, this study introduced a relief supply chain model under uncertainty with a leader-member exchange approach. Accordingly, on the upper level, the government, as the supplier, would provide distributors (governmental and nongovernmental) with relief goods so that they could distribute them to victims at demand points. Afterward, the model was solved using MOEA/D. According to the results, the selected algorithm has high accuracy and acceptable performance. Further, during each iteration, the objective functions became more optimized. It was also revealed that the objective function associated with the minimization of the unresponded demands outweighed the other objectives. Finally, the sensitivity of the objective functions to several parameters was examined. In the future research, we intend to locate distribution organizations based on appropriate criteria and also apply government budget restrictions and use simulation methods to extract better results.



---

*Source: 3352255-2023-02-20.xml*