# A Hybrid Mathematical and Decision-Making Model to Determine the Amount of Economic Order considering the Discount

Authors: Hanif Hazrati; Abbas Barzegarinegad; Hamid Siaby-Serajehlo
Journal: Mathematical Problems in Engineering
 (2021)

---

## Abstract

Suppliers are one of the most important parts of the supply chain, whose performance indirectly has a significant impact on customer satisfaction. Because customer demands are different from organizations, organizations have to consider different criteria for selecting their suppliers. In recent years, many studies in this field have been conducted using various criteria and methods. The main purpose defined in this research is to develop a model for simultaneous item ordering systems in real business conditions. In this research, a model is developed by considering the two objectives of minimizing overall costs and maximizing the amount of products ordered from different suppliers based on their weight value. Weights are calculated based on different criteria using the fuzzy analytic hierarchy process method for each supplier in different periods. Then, due to the multiobjective nature of the model, the proposed model has been solved by using the epsilon constraint in GAMS and nondominated sorting genetic algorithm II in MATLAB software. Considering the simultaneous order of inventory of multiproduct with several suppliers in several periods of time in discrete space with discount is one of the contributions of this research. To validate the proposed model, the results of the exact solution are compared with the meta-heuristic solution. Comparison results and assessment metrics indicate that the results of the proposed solution approach with an error of less than 1% had good performance. The results show that the system cost increases, by increasing the amount of discount, because of the increase in the amount of demand. Therefore, with a 30% increase in the discount, the system costs will increase to 36,496 units. Also, with a 20% reduction, the cost reduction will be reduced to 14,170 units.


---

## Body

## 1. Introduction

Choosing a set of suppliers is critical to the success of organizations. In recent years, more attention has been paid to the choice of suppliers [1, 2]. Effective selection and evaluation of suppliers is an important responsibility that should always be considered by purchasing managers [3, 4]. The criticality of supplier selection is due to the effects it has on the elements of organizations’ end products [5, 6]. According to Mahmoudi et al. [7], the decision to select a supplier involves the following three major decisions: (1) what products to order? Most papers written in this field are for situations where suppliers are selected for only one product. In these cases, different internal dependencies such as reductions in inspection and ordering costs that may exist for a variety of products are not considered. (2) In what quantity and from which supplier should the order be made? (3) In what time periods should the order be made? Determining the optimum quantity of inventory and selecting suppliers are closely related [8]. Making a decision on scheduling orders over time by choosing a supplier can dramatically reduce costs over time [9]. Therefore, one of the most important decisions is to buy, select, and maintain a complete group of suppliers, but recently another concept called supply chain management has emerged [10]. This has led researchers to increasingly realize that the issue of supplier selection is a factor that increases competition throughout the supply chain.In such issues, the fixed costs of ordering to the supplier and the bid price of that supplier have been considered. On the other hand, the issue of inventory management and control has been studied extensively in the literature [11, 12]. The issue of determining the amount of order is one of the important issues in this area. In these issues, two types of inventory costs and ordering costs have been considered [13]. By reducing the amount of goods in each order and as a result of increasing the number of orders, the inventory of goods in stock and related costs will decrease, but ordering costs will increase, which is not actually cost-effective. On the other hand, the cost of arranging an order from a supplier for a number of different products consists of two components: (1) the total cost of ordering is independent of the number of different products in an order and (2) the ordering cost that depends on the number of different products in an order [14]. The above problem is known as the problem of simultaneous completion of inventory. Due to the overall cost of ordering, using group ordering can lead to significant cost savings [15]. This savings increases significantly when demand between items is closely related [16]. There are two types of strategies for solving inventory completion problems: (1) direct grouping strategy and (2) indirect grouping strategy. In the direct grouping strategy, products are divided into a number of predefined groups and products within each group are ordered with a common cycle time. In the indirect grouping strategy, ordering is performed at regular intervals, which is usually called the base cycle time. Each product has an integer coefficient from the base cycle time and is ordered in custom cycles that is an integer multiple of this coefficient. Research questions include the following: (1) What providers are selected to serve? and (2) How can supply chain supply costs be reduced as much as possible?


## 2. Literature Review

Most of the decision models proposed for the supplier selection process have focused on the final supplier selection phase. Some of this research has been reviewed in this area. Gaballa [17] was the first researcher to apply mathematical programming to supplier selection in a real case. He has used an integer mixed programming model to formulate this decision in the Australian Post Office. Khalili Nasr et al. [18] proposed a multiobjective fuzzy model to minimize costs in the supply chain. The proposed two-echelon model selects suppliers and assigns them to manufacturers. The main purposes of this study are to minimize environmental costs, operating costs, and lost demand, and maximize employment. Finally, the proposed model is solved by goal programming. Li et al. [19] designed a green supply chain with a single commodity. The most important goals of this research include maximizing the profits of manufacturers, distributors, and collection centers. Another goal is to minimize environmental costs. Considering priority for customers is one of the contributions of this research. The results indicate that with increasing demand, environmental and transportation costs increase sharply. Diabat and Jebali [20] presented a multiperiod, multicommodity model in the supply chain. The intended products have different qualities, and based on this quality, recovery is done. There is also a special penalty for each return. The objectives are to minimize chain costs while minimizing environmental impact. The results of sensitivity analysis indicate that with increasing demand, the amount of shipping costs increases sharply. RezaHoseini et al. [10] have presented a linear programming model with several criteria. In this model, two sets of factors are considered: (1) supplier characteristics including quality, price, and delivery date and (2) characteristics of the buyer company including material requirements and assurance inventory. Shahrokhi et al. [9] proposed an ideal nonlinear mixed nonlinear model for supplier selection. Price, quality, delivery, and service are included in this model, and all criteria are considered in the presented format. Beiki et al. [21] developed a nonlinear mathematical modeling and a heuristic method for selecting suppliers under the terms of multiple goods, multiple suppliers, resource constraints, and discounts based on the number of purchases. Warehouse and investment constraint are considered a limitation in the model. Firouzi and Jadidi [22] used linear and complex number integer programming to select suppliers. In their model, price, delivery, quality, and discount are considered. The quality and time of delivery are considered a limitation. Babbar and Amin [23] used multi-objective linear programming to select suppliers. In this model, total price, quality, and late delivery are considered goals, and two sets of constraints are as follows: (1) system constraints (constraints that are not directly under the control of the purchasing manager) such as supplier capacity, minimum order quantity, and total budget purchasing and (2) policy restrictions including the maximum/minimum amount of purchases from a minimum/maximum number of suppliers to be employed. Baraki and Kianfar [24] proposed a mixed integer model to solve this problem in which buyers need to buy different items from different sellers whose capacity, quality, and quantity of deliveries are limited and the products are categorized. Ahmadi and Amin [25] also proposed a model to examine supplier selection in the presence of several suppliers, multiple criteria, and discounted prices. They have also considered the effects of budget constraints, quality, and capacity of suppliers. In the case of multiproduct issues, the issue is the issue of the simultaneous ordering of products, which occurs in one of the following cases: (1) when multiple products are ordered from one supplier, (2) when multiple products use common transportation equipment, and (3) when a product is packaged in different quantities after mass production or batch. Bektur [26] presented an integrated methodology for the selection of sustainable suppliers and order allocation problem. In their article, first, suppliers are ranked by using integrated the fuzzy analytic hierarchy process. Then, the fuzzy multiobjective optimization (FMOO) model is proposed to find quantities to be purchased from suppliers. The results indicate the proper performance of the proposed model. Kumar and Mahapatra [27] presented a multiwarehouse inventory model for an optimal replenishment policy. They solved their model by the rain optimization algorithm (ROA) considering different parameters such as total cost, total delivery time, and total investment on each item. The results indicate the proper performance of the proposed model. The contributions considered in this research are as follows: (1) considering the discount level in the proposed mathematical model, (2) providing a multiperiod and multicommodity model to control inventory and check the flow of products, and (3) providing a combined approach of nondominated sorting genetic algorithm II (NSGA-II) and fuzzy analytic hierarchy process (FAHP) to select suppliers in the proposed model.


## 3. Mathematical Modeling

The mathematical model presented in this article is an integer mathematical programming model. The purpose of presenting this model is to select suppliers dynamically in several discrete periods by considering the simultaneous ordering of items and discounts, which causes different suppliers to be selected in each period according to the limitations. Using the simultaneous ordering of items reduces the overall costs of sending products from different suppliers in different periods, and on the other hand, in the competitive market, each supplier tries to apply a discount in exchange for purchasing more of its competitive power. In fact, they have to balance the lower price benefits with more storage and transportation costs to make purchasing decisions, which makes it harder to decide. The assumptions of this issue include the following:(1)
Model parameters (purchase prices, demand, capacity, etc.) in different time periods are definite and known(2)
Shortage is not allowed in any time period(3)
Inventory holding cost has a linear relationship with inventory and is allocated only to the amount of inventory remaining at the end of each period(4)
The amount of inventory at the beginning of the first period and the end of the last period is zero(5)
There are several items of goods, suppliers, and time periods(6)
The buyer can buy the required amount of demand from each product in any time period from one or more suppliers(1)Min=∑i=1P∑j=1S∑k=1K∑t=1TXijktPrijkt+∑i=1P∑t=1ThitIit+∑j=1S∑t=1TAjtzjt+∑i=1P∑j=1S∑k=1K∑t=1TaijtYijkt,(2)Max=∑i=1P∑j=1S∑k=1K∑t=1TXijktWjt,(3)St.Iit=Iit−1+∑k=1K∑j=1SXijkt−Dit,∀i,t,(4)∑i=1P∑k=1KYijkt≤MZjt,∀j,t,(5)∑i=1PViIit−1+∑k=1K∑j=1SXijkt≤C,∀t,(6)∑i=1P∑j=1S∑k=1KXijktPrijkt≤Bt,∀t,(7)qijktYijkt≤Xijkt<qijk+1tYijkt,∀i,k,j,t,(8)∑k=1KYijkt≤1,∀i,j,t,(9)Ii0=0,∀i,(10)Iit=0,∀i,(11)Xijkt&Iit≥0,∀i,j,k,t,(12)Zjt&Yijkt=0or1,∀i,j,k,t.Equation (1), the objective function of the model, is to minimize the total costs over the entire time horizon, which consists of 4 components: the first part is related to the total cost of purchasing goods; the second part includes the cost of maintaining inventory at the end of period t; the third part is related to the total cost of ordering from supplier j in period t, which does not depend on the number of goods; and the fourth part deals with the costs of ordering, which depend on the number of goods.Constraint (2), the second function of the model, represents the sum of goods ordered from suppliers in different time periods based on the weight value of suppliers in different periods. Constraint (3) represents the inventory control between periods. Constraint (4) guarantees the selection of supplier j in time period t. Constraint (5) is related to the maximum storage capacity. Constraint (6) is related to the purchase budget in each time period. Constraint (7) guarantees that the quantity of the ordered goods will be at the discount level offered by the supplier. It is important to note that the upper limit of the last discount level offered by each supplier is equal to the production capacity of each supplier for the product i. Constraint (8) ensures that when ordering goods i, only one discount level is selected from the levels offered by supplier j. Constraints (9) and (10) show that inventory at the beginning and end of the time period for each commodity is zero.


## 4. Methodology

As stated in this study, a multi-objective model of linear programming integrated integer for supplier selection at different times according to the simultaneous ordering of products is presented. In this section, according to several objectives, considering the objective function, the NSGA-II meta-heuristic algorithm in the MATLAB software space and the epsilon-constraint method in the GAMS software space have been used to specify the part points of the model. The second objective function of the proposed model represents the value of the goods ordered from suppliers based on the criteria of quality, delivery of goods, and flexibility of service level, the weighted value for each supplier in various periods has been calculated using the fuzzy AHP method. In order to show the accuracy of the NSGA-II algorithm presented in this research, a number of problems in different dimensions have been solved using MATLAB software and compared with the results of the model using the epsilon-constraint method in GAMS.


### 4.1. Epsilon-Constraint Method (ECM)

The epsilon-constraint method is one of the well-known methods for dealing with multiobjective problems, which solves this type of problem by transferring all but one of the objective functions to the constraint at each stage. The Pareto boundary can be created by the epsilon-constraint method [12]:(13)minf1xs.t.x∈X,f2x≤ε2,fnx≤εn.The steps of the epsilon-constraint method are as follows:(1)
Select one of the objective functions as the main objective function(2)
Solve the problem according to one of the objective functions each time and get the optimal values of each objective function(3)
Divide the interval between the two optimal values of the subobjective functions into a predetermined number and obtain a table of values forε2,…,εn(4)
Each time obtain the problem with the main objective function with any of the valuesε2,…,εn(5)
Report the Parthian answers obtained


### 4.2. NSGA-II Algorithm

The NSGA-II method is a common method for solving problems with several objective functions based on genetic algorithms. This algorithm is an efficient method for solving problems with several objective functions, but it also has weaknesses for selecting dominant Pareto and in computational complexity. Therefore, a modified method called NSGA-II has been developed [28]. This method works better than the other methods because it uses the information SpandNp, the sum of the members of the population defeated by the P Pareto, and the number of times the P Pareto is defeated by other Pareto. In addition to all the functionality that NSGA-II has, it can be considered a model for the formation of many multi-objective optimization algorithms. This algorithm and its unique approach to multi-objective optimization problems have been used repeatedly by different people to create newer multi-objective optimization algorithms [29]. This algorithm solves the problem as it is defined and its purpose is not to turn the problem into a single objective problem. How to solve a problem with this algorithm requires a number of elements, one of which is the archive or archive of valid answers [30]. Because computer memory is limited, the number of archives cannot be allowed to grow as large as it would like. As a result, an archive size control mechanism is needed to control the number of archived responses.


#### 4.2.1. NSGA-II Algorithm Steps

To implement the NSGA-II algorithm, perform the following steps [31]:(1)
Generate a random initial solution of sizei = 1, ...., Popsize and set the number of iterations of the NSGA-II algorithm to 1.(2)
At this stage, the Pareto should be arranged on the basis of being defeated and they should be divided into fronts. The lower the number of fronts, the more Pareto the Pareto in them have defeated. To do this, combine the following steps for each of the Pareto, such asP:(2.1)
Considersp the sum of the members of a population defeated by the Pareto P, and consider its value to be zero.(2.2)
ConsiderNp the number of times the P Pareto is defeated by other Pareto and its value is zero.(2.3)
For each member of the populationn = 1, ..., Popsize as in q. Follow the steps below:(2.3.1)
If the Pareto P can defeat the Paretoq, then add q to the set sp.(2.3.2)
If theq Pareto can defeat the p Pareto, then add a unit to the set Np.(3)
If after examining all the Pareto points,Np=0, then it can be concluded that P is not defeated by any other Pareto points so P is added to the front F1.(4)
Continue all the following steps until the number of Pareto in fronti is equal to zero:(4.1)
Set the set of Pareto considered on thei + 1 front to Q and set it to zero, and then, perform the following steps for each P Pareto in the set fi:(4.2)
For each Pareto such asq in the set Sp on the fi front, perform the following steps:(4.2.1)
Subtract one unit fromNp. This indicates how many times the q Pareto has been defeated.(4.2.2)
IfNp equals zero, it indicates that the Pareto q is on the fi+1 front. In this case, you have to replace Q with q.(4.3)
Add a unit toi.(5)
After Pareto fronting, a number of them are selected based on the degree of defeat of other Pareto in order to create the next generation. The binary method is used to determine the Pareto. For this purpose, at the beginning, two Pareto points are randomly selected and a comparison is made between them, and whichever is better is added to the repository of answers. The criterion for better answers is based on the following two criteria:(5.1)
Rank priority: in this case, the answers that have a lower rank or a lower front are selected because the Pareto points of this front can dominate more Pareto points.(5.2)
In some cases, the two selected Pareto points may be in the same rank; in other words, they may both be in the same order. In this case, a criterion called CD is used, which is described below.(5.2.1)
For each frontni, consider fi the number of Pareto points in that front.(5.2.2)
Call the distance between the Pareto points on the frontsfi, and set the distance between all the Pareto points to zero.(5.2.3)
For each Pareto points such asj on the fi front, consider each of the objective functions of the problem as m and perform the following steps:(5.2.3.1)
On thefi front, arrange all the Pareto points based on the objective function m. In other words, we arrange the Pareto points on the fi front separately based on their objective functions.(5.2.3.2)
After arranging the Pareto points on thefi front based on the objective function m, set the distance between the first and last Pareto points to infinity. This is because there are no other Pareto points next to the Pareto points to cover it. For Pareto points 2 to n-1, the distance is obtained based on the following equations:(14)CDk=Idk1+…+Idkm,Idkm=Ik+1m−Ik−1mfmmax−fmmin.In the above equation, the congestion distance isIdkm in the objective function m. For the total congestion distance, Idk must be calculated and added for all objective functions, as specified in the above equation.(5.2.4)
After calculating the congestion distance, the Pareto points with the greater congestion distance are selected.(6)
After selecting the Pareto points, in the previous stage, a pond is created, which is called the selected population. Then, genetic operators are used to create the population of the children. The genetic operators used in this study are crossover and mutant operators.(7)
After determining the offspring ofPt, genetic operators, this population merges with the main population of Qt. Each pool has a capacity of n, and a number of Pareto points that have merged must be removed, so to reach capacity n, we must do the following steps:(7.1)
First, arrange the Pareto points according to the method described in step 2.(7.2)
Determine the distance of each Pareto points on the fronts.(7.3)
Start from frontf1 and select its Pareto points according to CD and pour into the new population pool (K + 1). Continue this step until the capacity of the new population pool (K + 1) reaches n. Figure 1 indicates the method of selecting a new population, crossover, and mutation.(8)
After creating the population (K + 1), go to step 2 and repeat these steps to the specified size.Figure 1
The proposed NSGA-II framework.


### 4.3. Calculate the Weight Importance of Each Supplier

As mentioned, the second objective function of the model includes maximizing purchases from suppliers based on the weight importance of suppliers in different time periods, which the weighted importance for each supplier is calculated using the fuzzy AHP method. Therefore, in this section, the method of calculating the weight of suppliers in different time periods is examined by the fuzzy AHP method.


#### 4.3.1. Fuzzy Analytical Hierarchical Process (FAHP)

Although the purpose of using the hierarchical analysis method is to obtain the opinion of experts and specialists, the analytical hierarchical process (AHP) method does not accurately reflect the way of human thinking because in this method, pairwise comparisons are considered numerically. In the fuzzy analytical hierarchical process (FAHP), after preparing the hierarchical diagram, decision-makers are asked to compare the elements of each level and express the relative importance of the elements using fuzzy numbers. The steps of the fuzzy analytical hierarchical process are as follows [32]:(1)
Drawing a hierarchical diagram(2)
Defining fuzzy numbers in order to perform even comparisons(3)
Forming an even comparison matrix using fuzzy numbersThe pairwise comparison matrix will be as follows:(15)A˜=1⋯a1n˜⋮⋱⋮an1˜⋯1.This matrix contains the following fuzzy numbers:(16)aij˜=1,i=j,1˜,3˜,5˜,7˜,9˜or1−1˜,3−1˜,5−1˜,7−1˜,9−1˜,i≠j.If the decision-making committee has several decision-makers, the components of the comprehensive pairwise comparison matrix used in the fuzzy hierarchical analytical method are a triangular fuzzy number whose first component is the minimum polls, whose second component is the average of the polls, and whose third component is maximum polls.(4)
CalculateSi for each row of the pairwise comparison matrix, where Si is a triangular fuzzy number and is calculated from the following equation: (17)Si=∑j=1mMgjj∗∑i=1n∑j=1mMgjj−1.In this regard,i represents the row number and j represents the column number. In this equation Mgjj is triangular fuzzy number. Also we can calculate ∑j=1mMgjj, ∑i=1n∑j=1mMgjj, and ∑i=1n∑j=1mMgjj−1 using the following equations:(18)∑j=1mMgjj=∑j=1mIj,∑j=1mmj,∑j=1muj,∑i=1n∑j=1mMgjj=∑i=1nIi,∑i=1mmi,∑i=1nui,∑i=1n∑j=1mMgjj−1=1∑i=1nui,1∑i=1nmi,1∑i=1nIi.In the above relations,Ii,mi,ui are the first to third components of fuzzy numbers.(5)
Calculate theSi relative to each other as follows:(19)VM2≥M1=1,ifm2≥m1,0,ifI1≥u2,I1−u2m2−u2−m1−I1,otherwise.(6)
Calculation of the weight of criteria and alternative in even comparison matrices is given as follows:(20)d′Aj=MinVSi≥Sk.Therefore, the non-normalized vector is as follows:(21)W′=d′A1,d′A2,…,d′AnT.(7)
Calculating the final weight vector as follows:(22)W=dA1,dA2,…,dAnT.


## 4.1. Epsilon-Constraint Method (ECM)

The epsilon-constraint method is one of the well-known methods for dealing with multiobjective problems, which solves this type of problem by transferring all but one of the objective functions to the constraint at each stage. The Pareto boundary can be created by the epsilon-constraint method [12]:(13)minf1xs.t.x∈X,f2x≤ε2,fnx≤εn.The steps of the epsilon-constraint method are as follows:(1)
Select one of the objective functions as the main objective function(2)
Solve the problem according to one of the objective functions each time and get the optimal values of each objective function(3)
Divide the interval between the two optimal values of the subobjective functions into a predetermined number and obtain a table of values forε2,…,εn(4)
Each time obtain the problem with the main objective function with any of the valuesε2,…,εn(5)
Report the Parthian answers obtained


## 4.2. NSGA-II Algorithm

The NSGA-II method is a common method for solving problems with several objective functions based on genetic algorithms. This algorithm is an efficient method for solving problems with several objective functions, but it also has weaknesses for selecting dominant Pareto and in computational complexity. Therefore, a modified method called NSGA-II has been developed [28]. This method works better than the other methods because it uses the information SpandNp, the sum of the members of the population defeated by the P Pareto, and the number of times the P Pareto is defeated by other Pareto. In addition to all the functionality that NSGA-II has, it can be considered a model for the formation of many multi-objective optimization algorithms. This algorithm and its unique approach to multi-objective optimization problems have been used repeatedly by different people to create newer multi-objective optimization algorithms [29]. This algorithm solves the problem as it is defined and its purpose is not to turn the problem into a single objective problem. How to solve a problem with this algorithm requires a number of elements, one of which is the archive or archive of valid answers [30]. Because computer memory is limited, the number of archives cannot be allowed to grow as large as it would like. As a result, an archive size control mechanism is needed to control the number of archived responses.


### 4.2.1. NSGA-II Algorithm Steps

To implement the NSGA-II algorithm, perform the following steps [31]:(1)
Generate a random initial solution of sizei = 1, ...., Popsize and set the number of iterations of the NSGA-II algorithm to 1.(2)
At this stage, the Pareto should be arranged on the basis of being defeated and they should be divided into fronts. The lower the number of fronts, the more Pareto the Pareto in them have defeated. To do this, combine the following steps for each of the Pareto, such asP:(2.1)
Considersp the sum of the members of a population defeated by the Pareto P, and consider its value to be zero.(2.2)
ConsiderNp the number of times the P Pareto is defeated by other Pareto and its value is zero.(2.3)
For each member of the populationn = 1, ..., Popsize as in q. Follow the steps below:(2.3.1)
If the Pareto P can defeat the Paretoq, then add q to the set sp.(2.3.2)
If theq Pareto can defeat the p Pareto, then add a unit to the set Np.(3)
If after examining all the Pareto points,Np=0, then it can be concluded that P is not defeated by any other Pareto points so P is added to the front F1.(4)
Continue all the following steps until the number of Pareto in fronti is equal to zero:(4.1)
Set the set of Pareto considered on thei + 1 front to Q and set it to zero, and then, perform the following steps for each P Pareto in the set fi:(4.2)
For each Pareto such asq in the set Sp on the fi front, perform the following steps:(4.2.1)
Subtract one unit fromNp. This indicates how many times the q Pareto has been defeated.(4.2.2)
IfNp equals zero, it indicates that the Pareto q is on the fi+1 front. In this case, you have to replace Q with q.(4.3)
Add a unit toi.(5)
After Pareto fronting, a number of them are selected based on the degree of defeat of other Pareto in order to create the next generation. The binary method is used to determine the Pareto. For this purpose, at the beginning, two Pareto points are randomly selected and a comparison is made between them, and whichever is better is added to the repository of answers. The criterion for better answers is based on the following two criteria:(5.1)
Rank priority: in this case, the answers that have a lower rank or a lower front are selected because the Pareto points of this front can dominate more Pareto points.(5.2)
In some cases, the two selected Pareto points may be in the same rank; in other words, they may both be in the same order. In this case, a criterion called CD is used, which is described below.(5.2.1)
For each frontni, consider fi the number of Pareto points in that front.(5.2.2)
Call the distance between the Pareto points on the frontsfi, and set the distance between all the Pareto points to zero.(5.2.3)
For each Pareto points such asj on the fi front, consider each of the objective functions of the problem as m and perform the following steps:(5.2.3.1)
On thefi front, arrange all the Pareto points based on the objective function m. In other words, we arrange the Pareto points on the fi front separately based on their objective functions.(5.2.3.2)
After arranging the Pareto points on thefi front based on the objective function m, set the distance between the first and last Pareto points to infinity. This is because there are no other Pareto points next to the Pareto points to cover it. For Pareto points 2 to n-1, the distance is obtained based on the following equations:(14)CDk=Idk1+…+Idkm,Idkm=Ik+1m−Ik−1mfmmax−fmmin.In the above equation, the congestion distance isIdkm in the objective function m. For the total congestion distance, Idk must be calculated and added for all objective functions, as specified in the above equation.(5.2.4)
After calculating the congestion distance, the Pareto points with the greater congestion distance are selected.(6)
After selecting the Pareto points, in the previous stage, a pond is created, which is called the selected population. Then, genetic operators are used to create the population of the children. The genetic operators used in this study are crossover and mutant operators.(7)
After determining the offspring ofPt, genetic operators, this population merges with the main population of Qt. Each pool has a capacity of n, and a number of Pareto points that have merged must be removed, so to reach capacity n, we must do the following steps:(7.1)
First, arrange the Pareto points according to the method described in step 2.(7.2)
Determine the distance of each Pareto points on the fronts.(7.3)
Start from frontf1 and select its Pareto points according to CD and pour into the new population pool (K + 1). Continue this step until the capacity of the new population pool (K + 1) reaches n. Figure 1 indicates the method of selecting a new population, crossover, and mutation.(8)
After creating the population (K + 1), go to step 2 and repeat these steps to the specified size.Figure 1
The proposed NSGA-II framework.


## 4.2.1. NSGA-II Algorithm Steps

To implement the NSGA-II algorithm, perform the following steps [31]:(1)
Generate a random initial solution of sizei = 1, ...., Popsize and set the number of iterations of the NSGA-II algorithm to 1.(2)
At this stage, the Pareto should be arranged on the basis of being defeated and they should be divided into fronts. The lower the number of fronts, the more Pareto the Pareto in them have defeated. To do this, combine the following steps for each of the Pareto, such asP:(2.1)
Considersp the sum of the members of a population defeated by the Pareto P, and consider its value to be zero.(2.2)
ConsiderNp the number of times the P Pareto is defeated by other Pareto and its value is zero.(2.3)
For each member of the populationn = 1, ..., Popsize as in q. Follow the steps below:(2.3.1)
If the Pareto P can defeat the Paretoq, then add q to the set sp.(2.3.2)
If theq Pareto can defeat the p Pareto, then add a unit to the set Np.(3)
If after examining all the Pareto points,Np=0, then it can be concluded that P is not defeated by any other Pareto points so P is added to the front F1.(4)
Continue all the following steps until the number of Pareto in fronti is equal to zero:(4.1)
Set the set of Pareto considered on thei + 1 front to Q and set it to zero, and then, perform the following steps for each P Pareto in the set fi:(4.2)
For each Pareto such asq in the set Sp on the fi front, perform the following steps:(4.2.1)
Subtract one unit fromNp. This indicates how many times the q Pareto has been defeated.(4.2.2)
IfNp equals zero, it indicates that the Pareto q is on the fi+1 front. In this case, you have to replace Q with q.(4.3)
Add a unit toi.(5)
After Pareto fronting, a number of them are selected based on the degree of defeat of other Pareto in order to create the next generation. The binary method is used to determine the Pareto. For this purpose, at the beginning, two Pareto points are randomly selected and a comparison is made between them, and whichever is better is added to the repository of answers. The criterion for better answers is based on the following two criteria:(5.1)
Rank priority: in this case, the answers that have a lower rank or a lower front are selected because the Pareto points of this front can dominate more Pareto points.(5.2)
In some cases, the two selected Pareto points may be in the same rank; in other words, they may both be in the same order. In this case, a criterion called CD is used, which is described below.(5.2.1)
For each frontni, consider fi the number of Pareto points in that front.(5.2.2)
Call the distance between the Pareto points on the frontsfi, and set the distance between all the Pareto points to zero.(5.2.3)
For each Pareto points such asj on the fi front, consider each of the objective functions of the problem as m and perform the following steps:(5.2.3.1)
On thefi front, arrange all the Pareto points based on the objective function m. In other words, we arrange the Pareto points on the fi front separately based on their objective functions.(5.2.3.2)
After arranging the Pareto points on thefi front based on the objective function m, set the distance between the first and last Pareto points to infinity. This is because there are no other Pareto points next to the Pareto points to cover it. For Pareto points 2 to n-1, the distance is obtained based on the following equations:(14)CDk=Idk1+…+Idkm,Idkm=Ik+1m−Ik−1mfmmax−fmmin.In the above equation, the congestion distance isIdkm in the objective function m. For the total congestion distance, Idk must be calculated and added for all objective functions, as specified in the above equation.(5.2.4)
After calculating the congestion distance, the Pareto points with the greater congestion distance are selected.(6)
After selecting the Pareto points, in the previous stage, a pond is created, which is called the selected population. Then, genetic operators are used to create the population of the children. The genetic operators used in this study are crossover and mutant operators.(7)
After determining the offspring ofPt, genetic operators, this population merges with the main population of Qt. Each pool has a capacity of n, and a number of Pareto points that have merged must be removed, so to reach capacity n, we must do the following steps:(7.1)
First, arrange the Pareto points according to the method described in step 2.(7.2)
Determine the distance of each Pareto points on the fronts.(7.3)
Start from frontf1 and select its Pareto points according to CD and pour into the new population pool (K + 1). Continue this step until the capacity of the new population pool (K + 1) reaches n. Figure 1 indicates the method of selecting a new population, crossover, and mutation.(8)
After creating the population (K + 1), go to step 2 and repeat these steps to the specified size.Figure 1
The proposed NSGA-II framework.


## 4.3. Calculate the Weight Importance of Each Supplier

As mentioned, the second objective function of the model includes maximizing purchases from suppliers based on the weight importance of suppliers in different time periods, which the weighted importance for each supplier is calculated using the fuzzy AHP method. Therefore, in this section, the method of calculating the weight of suppliers in different time periods is examined by the fuzzy AHP method.


### 4.3.1. Fuzzy Analytical Hierarchical Process (FAHP)

Although the purpose of using the hierarchical analysis method is to obtain the opinion of experts and specialists, the analytical hierarchical process (AHP) method does not accurately reflect the way of human thinking because in this method, pairwise comparisons are considered numerically. In the fuzzy analytical hierarchical process (FAHP), after preparing the hierarchical diagram, decision-makers are asked to compare the elements of each level and express the relative importance of the elements using fuzzy numbers. The steps of the fuzzy analytical hierarchical process are as follows [32]:(1)
Drawing a hierarchical diagram(2)
Defining fuzzy numbers in order to perform even comparisons(3)
Forming an even comparison matrix using fuzzy numbersThe pairwise comparison matrix will be as follows:(15)A˜=1⋯a1n˜⋮⋱⋮an1˜⋯1.This matrix contains the following fuzzy numbers:(16)aij˜=1,i=j,1˜,3˜,5˜,7˜,9˜or1−1˜,3−1˜,5−1˜,7−1˜,9−1˜,i≠j.If the decision-making committee has several decision-makers, the components of the comprehensive pairwise comparison matrix used in the fuzzy hierarchical analytical method are a triangular fuzzy number whose first component is the minimum polls, whose second component is the average of the polls, and whose third component is maximum polls.(4)
CalculateSi for each row of the pairwise comparison matrix, where Si is a triangular fuzzy number and is calculated from the following equation: (17)Si=∑j=1mMgjj∗∑i=1n∑j=1mMgjj−1.In this regard,i represents the row number and j represents the column number. In this equation Mgjj is triangular fuzzy number. Also we can calculate ∑j=1mMgjj, ∑i=1n∑j=1mMgjj, and ∑i=1n∑j=1mMgjj−1 using the following equations:(18)∑j=1mMgjj=∑j=1mIj,∑j=1mmj,∑j=1muj,∑i=1n∑j=1mMgjj=∑i=1nIi,∑i=1mmi,∑i=1nui,∑i=1n∑j=1mMgjj−1=1∑i=1nui,1∑i=1nmi,1∑i=1nIi.In the above relations,Ii,mi,ui are the first to third components of fuzzy numbers.(5)
Calculate theSi relative to each other as follows:(19)VM2≥M1=1,ifm2≥m1,0,ifI1≥u2,I1−u2m2−u2−m1−I1,otherwise.(6)
Calculation of the weight of criteria and alternative in even comparison matrices is given as follows:(20)d′Aj=MinVSi≥Sk.Therefore, the non-normalized vector is as follows:(21)W′=d′A1,d′A2,…,d′AnT.(7)
Calculating the final weight vector as follows:(22)W=dA1,dA2,…,dAnT.


## 4.3.1. Fuzzy Analytical Hierarchical Process (FAHP)

Although the purpose of using the hierarchical analysis method is to obtain the opinion of experts and specialists, the analytical hierarchical process (AHP) method does not accurately reflect the way of human thinking because in this method, pairwise comparisons are considered numerically. In the fuzzy analytical hierarchical process (FAHP), after preparing the hierarchical diagram, decision-makers are asked to compare the elements of each level and express the relative importance of the elements using fuzzy numbers. The steps of the fuzzy analytical hierarchical process are as follows [32]:(1)
Drawing a hierarchical diagram(2)
Defining fuzzy numbers in order to perform even comparisons(3)
Forming an even comparison matrix using fuzzy numbersThe pairwise comparison matrix will be as follows:(15)A˜=1⋯a1n˜⋮⋱⋮an1˜⋯1.This matrix contains the following fuzzy numbers:(16)aij˜=1,i=j,1˜,3˜,5˜,7˜,9˜or1−1˜,3−1˜,5−1˜,7−1˜,9−1˜,i≠j.If the decision-making committee has several decision-makers, the components of the comprehensive pairwise comparison matrix used in the fuzzy hierarchical analytical method are a triangular fuzzy number whose first component is the minimum polls, whose second component is the average of the polls, and whose third component is maximum polls.(4)
CalculateSi for each row of the pairwise comparison matrix, where Si is a triangular fuzzy number and is calculated from the following equation: (17)Si=∑j=1mMgjj∗∑i=1n∑j=1mMgjj−1.In this regard,i represents the row number and j represents the column number. In this equation Mgjj is triangular fuzzy number. Also we can calculate ∑j=1mMgjj, ∑i=1n∑j=1mMgjj, and ∑i=1n∑j=1mMgjj−1 using the following equations:(18)∑j=1mMgjj=∑j=1mIj,∑j=1mmj,∑j=1muj,∑i=1n∑j=1mMgjj=∑i=1nIi,∑i=1mmi,∑i=1nui,∑i=1n∑j=1mMgjj−1=1∑i=1nui,1∑i=1nmi,1∑i=1nIi.In the above relations,Ii,mi,ui are the first to third components of fuzzy numbers.(5)
Calculate theSi relative to each other as follows:(19)VM2≥M1=1,ifm2≥m1,0,ifI1≥u2,I1−u2m2−u2−m1−I1,otherwise.(6)
Calculation of the weight of criteria and alternative in even comparison matrices is given as follows:(20)d′Aj=MinVSi≥Sk.Therefore, the non-normalized vector is as follows:(21)W′=d′A1,d′A2,…,d′AnT.(7)
Calculating the final weight vector as follows:(22)W=dA1,dA2,…,dAnT.


## 5. Computational Results

In this section, the computational results obtained from solving the model in small, medium, and large scales will be shown. To demonstrate the efficiency of the NSGA-II algorithm, five small-scale problems are generated and the results are compared with the results of theε−constraint method. Table 1 shows the range of considered parameters. The parameter values follow a uniform distribution function. For example, ordering costs can be between 100 and 600 units.Table 1
Data ranges of parameters used in the model
ParametersRange of valuesPrijkt∼U150000,300000Vi∼U100000,200000Bt∼U100000,250000qijkt∼U200000,800000Lijt∼U100000,300000aijt∼U50,100Hit∼U50,100Ajt∼U100,600Table2 shows the number of products, suppliers, discount levels, and periods. As it is known, with increasing the scale of the problem, the number of products, suppliers, quality levels, and periods increase.Table 2
Dimensions of problem instances and sizes.
Problem sizeSizeijkPeriodSmallSS12323SS23433SS34643SS45884SS56894MediumMS1111194MS21213105MS31414125MS41515126MS51716136LargeLS12025205LS22227306LS32529408LS42830508LS53035609Tables3 and 4 show the results of the FAHP approach. As it is known, Table 3 calculates the incompatibility rate in pairs. This rate is calculated to be 0.78. Also, the weight of the first objective function is calculated according to experts 0.488 and the weight of the second objective function is equal to 0.512.Table 3
Incompatibility rate in pairs.
VSc21≥ScijVSc22≥Sci—10.78—Table 4
Weight of objective functions.
FAHPOBJ1OBJ2Final weight0.4880.512Rank21Table5 shows a comparison of the results of the epsilon-constraint and NSGA-II approaches. As it turns out, the results of the first and second objective functions are reported by both methods. The results show that the calculated error rate for all problems was less than 1%. Solution time also increases exponentially in the epsilon-constraint method, whereas it increases at a much slower rate in the NSGA-II method. Therefore, the solution results of the proposed NSGA-II approach can be trusted.Table 5
The results of calculation for small, medium, and large size test problems.
NoEpsilon constraintNSGA-IIError %f1f2Time (s)f1f2Time (s)f1f2135473.21260.6235473.21260.6200235498.61340.72635498.61340.7600335532.01489.95235541.51484.470.021435546.91553.89335560.31549.990.030.5535551.71675.514235576.21671.4140.061656788.82466.0136456799.92462.6230.011756798.42583.9289756823.22579.1260.040.07856836.22700.7374256849.52695.3320.020.9956864.42839.3498156891.12832.0380.040.91056892.72898.6642956938.72984.1410.0811173963.26783.2—73967.26780.0—0.0050.04711274556.26932.1—74559.96928.2—0.0040.0561375125.58988.0—75128.68984.4—0.0040.0401477068.59538.7—77074.89532.5—0.0080.0641579862.010725.0—79868.210720.0—0.0070.046Figure2 shows the Pareto solutions resulting from model solving. The blue points show the exact solution results, and the red points show the Meta-heuristic solution results. As can be seen, the results of exact solution and meta-heuristic are very close to each other, so the results of the meta-heuristic approach can be trusted.Figure 2
Pareto solutions.Figure3 shows a comparison of the solution times of the two approaches. As it is known, the solution time of epsilon-constraint approach increases exponentially, whereas the solution time of the NSGA-II approach is reasonable.Figure 3
CPU time.Table6 shows the assessment metrics results for the NSGA-II approach and the epsilon constraint. As can be seen, the average mean ideal distance for the epsilon-constraint and NSGA-II approaches is 6.491 and 6.552, respectively. Also, the average spacing metric for the epsilon-constraint and NSGA-II approaches is 0.401 and 0.404, respectively.Table 6
Multi-objective metrics obtained for each algorithm.
NoMidSMEpsilon constraintNSGA-IIEpsilon constraintNSGA-II16.236.250.3650.36526.266.280.3690.37136.336.350.3710.37446.366.380.3750.37756.426.460.3800.38666.606.650.4240.42876.636.670.4270.43086.656.690.4300.43596.696.710.4360.438106.746.780.4390.442Figure4 shows the MID and SM metrics results. The blue points show the epsilon-constraint results, and the red points show the NSGA-II results. As can be seen, the results of exact solution and meta-heuristic are very close to each other, so the results of the meta-heuristic approach can be trusted.Figure 4
MID and SM metrics.Table7 shows the results the one of the Pareto points. As can be seen in Table 7 in period 1, suppliers 1, 3, 4, and 5 are active. Also in period 2, suppliers 1, 2, 3, 4, and 6 are active.Table 7
Supplier selection outputs.
Period 1Period 2Supplier 111Supplier 201Supplier 311Supplier 411Supplier 510Supplier 601Figure5 shows the effect of the change in discount on the cost of the proposed systems. As it has been determined, the system cost increases, by increasing the amount of discount, because of the increase in the amount of demand. Therefore, with a 30% increase in the discount, the system costs will increase to 36,496 units. Also, with a 20% reduction, the cost reduction will be reduced to 14,170 units.Figure 5
Effect of the change in discount on the cost.


## 6. Conclusion

One of the important factors for survival in today’s highly competitive environment is the reduction of production costs. Choosing the right suppliers can significantly reduce purchasing costs and increase an organization’s competitiveness. The reason for this is that in most industries, the cost of raw materials and components of the product accounts for a large part of the cost of the product. According to the obtained results, the average error rate between the proposed NSGA-II algorithm and the epsilon-constraint method is acceptable. It is less and more acceptable. The incompatibility rate in pairs is calculated to be 0.78 using FAHP. Also, the weight of the first objective function is calculated according to experts 0.488 and the weight of the second objective function is equal to 0.512. The average mean ideal distance for the epsilon-constraint and NSGA-II approaches is 6.491 and 6.552, respectively. Also, the average spacing metric for the epsilon-constraint and NSGA-II approaches is 0.401 and 0.404, respectively. The results show that the system cost increases, by increasing the amount of discount, because of the increase in the amount of demand. Therefore, with a 30% increase in the discount, the system costs will increase to 36,496 units. Also, with a 20% reduction, the cost reduction will be reduced to 14,170 units. Among the research limitations is the lack of access to accurate information on transportation costs. Therefore, transportation costs are based on the estimates of transportation specialists. To develop the proposed model, the following future suggestions can be considered:(1)
Considering the policy of shortage of products in different time periods(2)
Considering probability or fuzzy uncertainty for capacity and demand parameters(3)
Considering the time value of money in the issue of choosing a supplier and ordering items at the same time(4)
Considering the maintenance cost of each unit of products as decreasing or increasing over time (dynamic nature of maintenance cost)



---

*Source: 5229949-2021-12-02.xml*