A hybrid vendor managed inventory and redundancy allocation optimization problem in supply chain management: An NSGA-II with tuned parameters

Computers & Operations Research 41 (2014) 53–64

Contents lists available at ScienceDirect

Computers & Operations Research journal homepage: www.elsevier.com/locate/caor

A hybrid vendor managed inventory and redundancy allocation optimization problem in supply chain management: An NSGA-II with tuned parameters Javad Sadeghi a,n, Saeid Sadeghi b, Seyed Taghi Akhavan Niaki c a

Young Researchers and Elite Club, Qazvin Branch, Islamic Azad University, Qazvin, Iran Young Researchers and Elite Club, South Tehran Branch, Islamic Azad University, Tehran, Iran c Department of Industrial Engineering, Sharif University of Technology, Tehran, Iran b

art ic l e i nf o

a b s t r a c t

Available online 11 August 2013

In this research, a bi-objective vendor managed inventory model in a supply chain with one vendor (producer) and several retailers is developed, in which determination of the optimal numbers of different machines that work in series to produce a single item is considered. While the demand rates of the retailers are deterministic and known, the constraints are the total budget, required storage space, vendor's total replenishment frequencies, and average inventory. In addition to production and holding costs of the vendor along with the ordering and holding costs of the retailers, the transportation cost of delivering the item to the retailers is also considered in the total chain cost. The aim is to find the order size, the replenishment frequency of the retailers, the optimal traveling tour from the vendor to retailers, and the number of machines so as the total chain cost is minimized while the system reliability of producing the item is maximized. Since the developed model of the problem is NP-hard, the multiobjective meta-heuristic optimization algorithm of non-dominated sorting genetic algorithm-II (NSGA-II) is proposed to solve the problem. Besides, since no benchmark is available in the literature to verify and validate the results obtained, a non-dominated ranking genetic algorithm (NRGA) is suggested to solve the problem as well. The parameters of both algorithms are first calibrated using the Taguchi approach. Then, the performances of the two algorithms are compared in terms of some multi-objective performance measures. Moreover, a local searcher, named simulated annealing (SA), is used to improve NSGA-II. For further validation, the Pareto fronts are compared to lower and upper bounds obtained using a genetic algorithm employed to solve two single-objective problems separately. & 2013 Elsevier Ltd. All rights reserved.

Keywords: Vendor managed inventory Redundancy allocation NSGA-II NRGA Taguchi method

1. Introduction Supply chain management (SCM) plays an important role in competing firms of today's market. Such an integrated supply chain (SC) can reduce total cost compared to the cost when each part decides independently. The SCM is the coordination between location, inventory, transportation, and production for a set, which consists of a network of facilities and distribution options, in order to reach the best mix of efficiency and responsiveness for the market being served [1]. One of the most important issues in SCM is the inventory management, where the total inventory cost is the main cost and many researchers developed various inventory models and methodologies that lead into different managerial toolboxes. Strategies

n

Corresponding author. Tel.: +98 2155 912 155; fax: +98 2155 912 155. E-mail addresses: [email protected], [email protected] (J. Sadeghi), [email protected] (S. Sadeghi), [email protected] (S.T.A. Niaki). 0305-0548/$ - see front matter & 2013 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.cor.2013.07.024

such as quick response (QR), advanced continuous replenishment (ACR), and vendor-managed inventory (VMI) are a few examples of retailer–supplier partnerships that are to be considered in making an inventory-related decision. Among these strategies, VMI due to its better performances has recently become popular and many successful retailers like Wal-Mart, Kmart, and JcPenney used it because of its benefits [2,3]. In the VMI model, there is a cooperation between vendor and retailers to determine the inventory level. The vendor, who knows retailer's sales and inventory information, determines the replenishment frequency and the order size of the retailer. The most important decisions to be made in the VMI model are determining the order frequencies and the volumes of retailer's replenishments. From a general point of view, the VMI leads into a reduction in demand variability, and thus a decrease in the average inventory, which reduces the total inventory cost. In the VMI strategy, an approach that has been recently discussed to the ownership stock between the vendor and the retailer is the consignment stock (CS). In the CS, the vendor receives the money only when the retailer sells the goods. Although the retailer holds the

54

J. Sadeghi et al. / Computers & Operations Research 41 (2014) 53–64

stock, but the vendor is legally the owner of the unsold stock. An example of a work in this area is Zavanella and Zanoni [4] who proposed a VMI model for a supply chain consisting a single vendor and several retailers. Two issues seem worthy of investigation in a VMI model developed in the CS environment. They are


Since the vendor produces goods and delivers them to retailers,




it seems reasonable to find the optimal route that the vendor's vehicle must pass through it to deliver the goods to the retailers (similar to the traveling salesman problem (TSP)). In order to provide on-time deliveries of the goods produced by the vendor using several machines that work in series, it seems reasonable to determine the optimum number of each machine in parallel coordination such that the total production system reliability is maximized (similar to the redundancy allocation problem (RAP)). In the reliability theory, system reliability, which is the probability of a system successfully performs as designed, can be improved by redundancy. Redundancy is a design technique to enhance system reliability; in other words, redundancy increases system reliability with adding machines to the system regarding constraints such as available budget and space [5].

In this paper, the VMI model of Zavanella and Zanoni [4] is extended to include the above two scenarios. In other words, a biobjective VMI of a supply chain consisting of a single vendor and multiple retailers is proposed, in which the first objective is the minimization of the total chain cost including production, ordering, holding, and transportation. The second objective is the maximization of the production system reliability using the RAP approach. The aim is to find the order size, replenishment frequency of the retailers, the routing tour, and the number of machines of different types in series so as the above two objectives are simultaneously optimized. In addition, there are five constraints. First, the retailers' warehouse has a predetermined limited space; second, the average inventory level of the vendor is restricted to an upper bound; third, the replenishment frequency is limited; fourth, the space required for machine installation is constrained; and fifth there is an upper bound for the budget to privide the machines. Moreover, the vendor determines the shortest route to deliver goods so as to reduce the transportation cost. After showing that the developed model of the problem is NP-hard, a multi-objective meta-heuristic algorithm is developed in this research to find a near-optimum solution of the problem. Since there is no benchmark available in the literature to verify and validate the solution obtained, another multi-objective meta-heuristic algorithm is developed to find a near-optimum solution as well. Moreover, the Taguchi method adjusts the parameters of the algorithms to have solutions with better quality. In addition, a local searcher named simulated annealing (SA) improves meta-heuristic algorithms in order to make them faster. This research is most probably one of the pioneer works that investigates both the transportation cost and the redundancy allocation optimization for productive environments in the VMI models. The reminder of the paper is organized as follows. A literature review is brought in the next section. The mathematical formulation of the problem is introduced in Section 3. The solution algorithms and the parameter tuning of the algorithms are presented in Sections 4 and 5, respectively. The algorithms are evaluated and are improved in Section 6. Finally, conclusion and recommendations for future research are in Section 7.

economic order quantity (EOQ) model, the economic production quantity (EPQ) model, and solution algorithms proposed to solve VMI problems. These classes are discussed in the following three subsections, respectively. Moreover, since this research concerns with a VMI based on the EPQ, the majority of the research works are considered for the second class. 2.1. The VMI model based on the EOQ policy For the VMI models developed based on the EOQ policy, Yao et al. [6] were the first who proposed a mathematical formulation for onevendor and one-retailer SC showing the advantages of implementing the VMI approach. Then, Darwish and Odah [7] extended the model for one-vendor multi-retailers and provided a heuristic algorithm in which the computational effort was reduced. Next, Pasandideh et al. [8] suggested a genetic algorithm (GA) with calibrated parameters to find near-optimum solution of the one-supplier one-retailer SC problem with several constraints in which shortages were backordered. Afterward, Cárdenas et al. [9] improved the proposed VMI model of Pasandideh et al. [8] considering the available inventory in limited storage space environments. They [9] first presented a heuristic algorithm to solve the VMI problem, and then showed that their algorithm was better than the GA of Pasandideh et al. [8]. Although Cárdenas et al. [9] and Pasandideh et al. [8] studied the VMI policy under some specific assumptions, the replenishment frequency was not considered in their models. Therefore, Yao et al. [6] and Darwish and Odah [7] investigated a one-vendor oneretailer VMI model with replenishment frequency in which the vendor sends a lot of size Q to a retailer in n shipments of size q each. Moreover, Sadeghi et al. [10] proposed a parameter-tuned GA using response surface methodology (RSM) to find a near-optimum solution of the VMI problem under the above consideration. Hariga et al. [11] first extended the VMI model of Darwish and Odah [7] to include unequal replenishment intervals. Then, they used a simple heuristic procedure to solve the problem. Considering stochastic demands, Kiesmuller and Broekmeulen [12] presented a multiproduct VMI model to reduce the handling cost at a warehouse and the transportation cost between the warehouse and the stores. 2.2. The VMI model based on the EPQ policy In an EPQ policy where there is a production constraint, Goyal [13] presented an inventory model for the first time in which the vendor produces items in order to fulfill a retailer's order in lot sizes of an integer multiple of the order sizes. Then, He [14] extended the model for the multi-retailer case. However, the effect of the production rate on the average inventory was not considered in the modeling. Next, Lu [15] suggested a heuristic approach to find a near-optimum solution of a one-vendor multi-retailer VMI model based on the EPQ policy. While Goyal [13] and Lu [15] presented VMI models based on the EPQ policy in which the shipments are of a unique size, Hill [16] developed the model under different shipment sizes. Moreover, Hill and Omar [17] first reviewed several previous works based on the unit holding cost between the vendor and the buyer, and then investigated the problem under various shipment sizes. Regarding the consignment stock (CS) approach, Braglia and Zavanella [18] extended the model proposed by Hill [16] for a single-vendor single-buyer SC in production environments. More extension is due to Zavanella and Zanoni [4] who developed the model to the multi-buyer case, where a sensitivity analysis was presented for model parameters.

2. Literature review

2.3. Solution algorithms of the VMI models

The literature survey of this research on the vendor managed inventory (VMI) problem is divided into three classes based on the

Evolutionary approaches have become popular to solve various complex optimization problems in recent years. In supply-chain

J. Sadeghi et al. / Computers & Operations Research 41 (2014) 53–64

environments for example, some research works including [8,19–22] presented integer nonlinear programming models and proposed GA to find near optimal solutions of the problems. Besides, one of the important advantages of using meta-heuristics in finding feasible solution is to avoid being trapped in local optima. This benefit can also be used in methods used in numerical analysis. For example, Yu et al. [23] took advantage of a meta-heuristic, namely genetic algorithm (GA), for a Newton–Raphson method in order to optimize a VMI problem. In this problem, a vendor replenishes a deteriorating product to some retailers. Moreover, meta-heuristics can be mixed with other methods to optimize complex problems in combinational optimization. For instance, Yu et al. [24] combined a GA with dynamic programming algorithm to solve a retailer selection problem in a VMI system. In addition, Kuo and Han [25] considered a supply chain distribution problem and proposed GA and particle swarm optimization (PSO) as solution algorithms. Although evolutionary algorithms were applied for multi-objective optimization problems, there is just one research available in the literature [26] that presented a non-dominated sorting genetic algorithm (NSGA-II) with un-calibrated parameters as a multi-objective evolutionary algorithm in VMI environment.

where z A f1; 2; …; r g. Based on (1), the order quantity of the ith retailer is obtained by qi ¼

di n1 q1 d1 ni

r

Q v ¼ ∑ ni qi Regarding Eq. (2), Q V can be written as follows: r

Q v ¼ ∑ ðdi n1 q1 =d1 Þ

Fig. 1 presents the consumption of stock for a one-vendor and three-retailer SC with n1 ¼ 2, n2 ¼ 3 and n3 ¼ 2so that the vendor replenishes retailers seven times. With the above assumption, the total ordering cost (TOC v ) and the total holding cost (THC v ) of the vendor per time unit are given in (5) and (6), respectively:

3. Mathematical model The following notations and assumptions are used throughout this paper. an index used for a retailer; i ¼ 1; 2; :::; r an index used for a machine; j ¼ 1; 2; :::; m ordering costs for retailer i ordering cost of the vendor holding cost for retailer i holding cost of the vendor vendor's production rate (item/time unit) vendor replenishment frequency of retailer i (per production cycle time) (decision variable) vendor replenishment frequency of retailer 1 (per production cycle time) demand rate of retailer i (item/time unit) demand rate of the 1st retailer demand rate of the vendor order quantity for retailer i order quantity for retailer 1 (decision variable) total vendor's order quantity upper bound for average inventory of the vendor space required storing one unit of the product total available storage space for retailers upper bound for replenishment frequencies total cost of the VMI system production system reliability the reliability of machine j number of redundant machine j (decision variable) purchasing cost of machine j total available space to install machines occupied space by a machine available budget to install machines transportation cost from retailer k to retailer l production cycle time is 1 if retailer k is reached from retailer l, 0 otherwise (decision variable) Assuming an equal consumption interval for retailers, we have ni qi nz qz ¼ di dz

ð4Þ

i¼1

i¼1

di d1 D qi q1 Qv Z f F λ TC R Rj Nj Cj e s B Skl T xkl

ð3Þ

i¼1

r

n1

ð2Þ

where d1 is utilized to simplify modeling. Considering the replenishment frequency of retailers (ni ) by the vendor per production cycle time, total vendor's order quantity is

TOC v ¼ ∑ Adi =ðni qi Þ

i j Ai A hi H P ni

55

ð1Þ Fig. 1. The graphical representation of inventory levels [4].

ð5Þ

56

J. Sadeghi et al. / Computers & Operations Research 41 (2014) 53–64 r

THC v ¼ ∑ Hqi di =ð2PÞ

ð6Þ

i¼1

Moreover, the total ordering cost (TOC r ) and the total holding cost (THC r ) per time unit for the retailers are obtained by TOC r ¼ ∑ Adi =qi

ð7Þ




r




1 d ∑ h n q 1 i 2i¼1 i i i P

 þ

qi di P

 ð8Þ

ð9Þ

There are several constraints on this VMI system; first, the warehouse space for all retailers is limited to F. Second, vendor's average inventory restricts to Z. Third, replenishment frequency is limited to λ. Consequently, we have

  r d dq ð10Þ ∑ f ni qi 1 i þ i i r F P P i¼1 qi di rZ i ¼ 1 2P

ð11Þ

r

∑ ni r λ

r

8l ¼ 1 : r

ð21Þ

r

∑ xkl ¼ 1;

8k ¼ 1 : r

ð22Þ

Skl ¼ 1

for all k ¼ l

ð23Þ

In order to enhance the production system reliability, R, in which several machines work in series, the redundancy allocation problem is considered in this research. Thus, the second objective function along with its corresponding constraints that are provided to determine an optimal combination of the number of the machines, Nj , is i m h ð24Þ MaxR ¼ ∏ 1ð1Rj ÞNj j¼1

m

m

j¼1

j¼1

∑ C j Nj r B þ ∑ C j

ð25Þ

ð12Þ m

i¼1

In order to consider the transportation cost, tc, the following equation is added to TC V MI : r

tc ¼ ∑ ∑ Skl xkl

ð13Þ

k¼1l¼1

s.t. for all k ¼ l

ð14Þ

where Skl is the transportation cost from retailer k to retailer l. Similar to the traveling salesman problem (TSP), it is assumed that the vendor delivers the products to all the retailers using a single vehicle. It is also assumed that the optimal traveling tour from the vendor to retailers is a base to calculate transportation cost. In addition, the vehicle is not allowed to pass a retailer's location twice. Therefore, r

∑ xkl ¼ 1;

8l ¼ 1 : r

k¼1

∑ Nj s r e

ð26Þ

j¼1

where inequality (25) presents the budget constraint and inequality (26) illustrates that the occupied space by redundant machines is restricted to a predetermined value. Note that N j ; ni ; q1 ; andxkl all are non-negative integers. Interested readers are referred to handbook [27] for the redundancy allocation problem. Minimizing the total system inventory cost along with maximizing the total system reliability, the multi-objective VMI problem can then be modeled by 2 Ad1 Hn1 q1 r di d1 r MinTC ¼ þ þ ∑ ∑ Ai ni n1 q1 2d1 P i ¼ 1 ni n1 q1 i ¼ 1
 r n1 q1 d d þ ∑ h d 1 i þ i 2d1 i ¼ 1 i i P ni P r

r

m

þ ∑ ∑ Skl xkl þ ∑ C j N j ; j¼1

k¼1l¼1

ð15Þ

m

h

Nj

Max R ¼ ∏ 1ð1Rj Þ

i

j¼1

r

∑ xkl ¼ 1;

8k ¼ 1 : r

l¼1

ð16Þ

The purchasing cost of the production machines as a representative for the vendor's operating cost is another cost that should be considered in the objective function. It can be easily obtained by ∑m j ¼ 1 C j Nj . Based on what was assumed and derived above, the first objective function becomes r

2 di

r

r

r

m

k¼1l¼1

s.t. ∑ f di n1 q1





 r d d ∑ f di n1 q1 1 i þ i =d1 rF P ni P i¼1 ! r

n1 q1

2

∑ di =ni =2d1 P r Z

i¼1

r

∑ ni r λ

i¼1 m

m

j¼1

j¼1

m

∑ Nj s r e

j¼1 r

ð17Þ

j¼1

 d d 1 i þ i =d1 r F P ni P

s.t.

∑ C j Nj rB þ ∑ C j

Ad1 Hn1 q1 d1 þ þ ∑ An ∑ n1 q1 2d1 P i ¼ 1 ni n1 q1 i ¼ 1 i i
 n1 q1 r d d þ ∑ h d 1 i þ i 2d1 i ¼ 1 i i P ni P þ ∑ ∑ Skl xkl þ ∑ C j Nj ;

i¼1

ð20Þ

s.t.

r

∑

r

ð19Þ

r

l¼1

TC V MI ¼ TOC r þ THC r þ TOC v þ THC v

MinTC ¼

=2d1 P rZ

∑ ni r λ

k¼1

Besides, similar to Zavanella's model [4], the average total cost of the VMI system per time unit is

Skl ¼ 1

∑

i¼1

∑ xkl ¼ 1;

i¼1

r

n1 q1

2 di =ni

i¼1

r

THC r ¼

!

r

∑ xkl ¼ 1;

8l ¼ 1 : r

∑ xkl ¼ 1;

8k ¼ 1 : r

k¼1 r l¼1

ð18Þ

Skl ¼ 1 for all k ¼ l N j ; ni ; q1 ; xkl are non  negative integers

ð27Þ

J. Sadeghi et al. / Computers & Operations Research 41 (2014) 53–64

57

Since the problem modeled in (27) is of a bi-objective integer nonlinear programming type, in the next section, a meta-heuristic algorithm is presented to solve it.

4. A solution algorithm A main branch in the theory of computation, named computational complexity, considers classifying computational problems regarding their inherent difficulty. There are four classes (P, NP, NP-complete, and NP-hard) with respect to the difficulty level of solving problems. In this scope, there are several open questions such as Is P ¼NP problem or not? Besides, many real-world optimization problems belong to the class of NP-hard and in order to solve NP-hard problems, there are not provably efficient algorithms, i.e. exact methods cannot solve the problems in normal time. According to the performed studies, meta-heuristic algorithms are suitable tools to optimize this class of problems [28]. The Hamiltonian cycle problem, one of real-world optimization problems, was shown NP-complete by Richard M. Karp in 1972, which infers NP-hardness of TSP [29]. In addition, the redundancy allocation problems in reliability area are another example of NP-hard problems [30]. Chern [31] proved that a RAP of a simple system belongs to the class of NP-hard problems. Moreover, INLP problems are also NP-complete [32]. Since the proposed bi-objective model shown in (27) consists of the above three NP-hard problems, it is NP-hard as well. This justifies the use of a meta-heuristic algorithm. In this section, the well-known multiobjective evolutionary algorithm (MOEA) of NSGA-II is presented to solve the problem.

4.1. Non-dominated sorting genetic algorithm-II (NSGA-II) A population-based search MOEA can present a set of Pareto optimal solutions of multi-objective optimization problems involving two or more conflicting objectives. One of these MOEAs that is frequently used in many optimization problem as the best technique to generate Pareto frontiers is the non-dominated sorting genetic algorithm-II (NSGA-II) proposed by Deb et al. [33]. To start NSGA-II, one first randomly generates a population P 1 with size N p chromosomes (solutions) and then sorts the chromosomes in P 1 into several fronts of non-dominated solutions. Considering the obtained chromosomes using the tournament selection operator for P 1 , the offspring population O1 is created with respect to the crossover rate (P c ) and the mutation rate (P m ). Moreover, the algorithm obtains the cost and the coverage of each chromosome in P 1 and O1 . After merging P 1 and O1 to form Rt , the algorithm sorts Rt in several non-dominated fronts F i , where the best F i s form the next population P tþ1 . Since the size of P tþ1 is equal to the size of Pt, all of elements from a front cannot be in P tþ1 . Hence, when a front is added to P tþ1 incompletely, the crowding distance approach is applied. Consequently, the required population is organized from the top elements of the front without loosing good solutions (elitism). The algorithm creates Otþ1 from P tþ1 using a crowded tournament algorithm and crossover and mutation operators. Regarding the stopping criteria and iterating the above stages, the algorithm hopefully presents the best Pareto optimal solutions. In this paper, in order to stop the algorithm, a fixed number of generation determined by the Taguchi method presented in Section 5 is used. In addition, the death penalty approach is used to penalize infeasible solutions that do not satisfy the constraints. Figs. 2 and 3 show a graphical representation of NSGA-II. For more details on the implementation of NSGA-II refer to [33,34].

Fig. 2. Flowchart of NSGA-II.

4.2. Non-dominated ranking genetic algorithms (NRGA) Since no benchmark is available in the literature to verify and validate the solutions obtained by NSGA-II, in this section a second popular MOEA called non-dominated ranking genetic algorithm (NRGA) is developed to obtain Pareto fronts. In addition, the Pareto optimal fronts are compared to the single-objective solutions obtained by a genetic algorithm (GA) applied to two singleobjective problems, separately. The single-objective solutions can be treated as lower and upper bounds to compare the Pareto fronts. Jadaan et al. [35] presented NRGA by exchanging the selection strategy of NSGA-II from the tournament selection to the roulette wheel. In the selection strategy of NRGA, there are two rankedbased roulette wheel tires; the first tier is to select the front, and the second tier is to select solution from the front. Eq. (28) gives the probability of selecting the ith front: Pi ¼

2  Ranki ; N F ðNF þ 1Þ

8 i ¼ 1; …; N F

ð28Þ

Here N F is the number of fronts and Ranki is the rank of front i. It can be understood from Eq. (28) that in order to select, the highest rank has the highest chance (in case of maximization). In addition, the probability of selecting the jth individual of the ith front is P ij ¼

2  Rankij ; I i ðI i þ 1Þ

8 i ¼ 1; …; NF ;

8 j ¼ 1; …; I i

ð29Þ

58

J. Sadeghi et al. / Computers & Operations Research 41 (2014) 53–64

Fig. 3. Graphical representation of NSGA-II [34].

Table 4 Normalized metrics obtained by NSGA-II.

Fig. 4. Chromosome representation.

Table 1 NSGA-II and NRGA parameters. Variable

Value

Np It Pc Pm

55 400 0.6 0.1

60 700 0.7 0.15

65 1000 0.8 0.2

NPS

GD

Sp

Best_Sol

Sum

0.3171 0.7317 0.9268 0.5854 0.2683 0.7561 0.4146 1.0000 1.0000

0.3300 0.7122 0.8322 0.5423 0.2849 0.8379 0.4796 1.0000 0.8833

0.6378 0.5680 0.7705 0.5978 0.2766 0.7739 0.8882 1.0000 0.8266

0.8353 0.9236 0.9531 0.9419 0.8264 0.8882 0.8455 0.9261 1.0000

2.1202 2.9356 3.4826 2.6673 1.6561 3.2561 2.6280 3.9261 3.7100

Table 5 Experimental results of NSGA-II. Table 2 Basic data of the numerical example with four retailers and seven machines. Rj

Cj

di

Ai

hi

0.65

19342

8361

58

6

0.25 0.36 0.11 0.73 0.95 0.77

11293 33194 43688 32276 43596 18198

7376 1388 2313

95 27 68

8 2 7

Transportation costs between retailers

R1 R2 R3 R4

R1

R2

R3

R4

0 473 319 576

473 0 399 708

319 399 0 676

567 708 676 0

H¼ 12; P ¼23974; A ¼ 265; B¼ 450000; s¼ 3; e¼100; Z ¼2000; f¼ 3; F¼ 3000; λ ¼120.

Table 3 The metrics obtained by NSGA-II. NPS

GD

Sp

Best_Sol

13 30 38 24 11 31 17 41 41

32.9515 15.2664 13.0660 20.0517 38.1710 12.9768 22.6723 10.8732 12.3094

6174.6608 6933.3911 5111.1861 6588.1991 14240.0623 5088.9630 4433.9113 3938.3916 4764.4153

0.8444 0.7637 0.7401 0.7489 0.8536 0.7942 0.8343 0.7617 0.7054

where I i is the number of individuals in front i. NRGA ranks the fronts and the individuals in a front based on the non-dominated rank and their crowding distance, respectively. After calculating P i

Np

It

Pc

Pm

Sum1

Sum2

Sum3

Sum4

Sum5

S/N

1 1 1 2 2 2 3 3 3

1 2 3 1 2 3 1 2 3

1 2 3 2 3 1 3 1 2

1 2 3 3 1 2 2 3 1

2.1202 2.9356 3.4826 2.6673 1.6561 3.2561 2.6280 3.9261 3.7100

1.2230 2.9122 3.0602 1.5643 3.4900 3.2422 1.9715 3.2697 3.8857

2.1580 2.3389 3.4002 2.5411 1.8442 3.8265 2.4219 3.4822 3.1287

1.9548 2.0107 2.7098 1.7041 3.3349 3.9820 2.3518 3.5118 3.4272

2.0467 2.1858 3.3078 3.0151 2.9452 2.9463 2.2918 3.6071 3.2671

4.9319 7.5775 9.9714 6.3510 7.2310 10.5947 7.2411 10.9819 10.7586

Table 6 Experimental results of NRGA. Np

It

Pc

Pm

Sum1

Sum2

Sum3

Sum4

Sum5

S/N

1 1 1 2 2 2 3 3 3

1 2 3 1 2 3 1 2 3

1 2 3 2 3 1 3 1 2

1 2 3 3 1 2 2 3 1

2.3784 1.8354 3.6135 2.0599 2.4176 3.4687 2.6394 3.4356 2.8538

2.1264 1.4594 2.5086 1.9432 1.8162 2.6866 1.9124 3.1035 1.7903

1.5834 2.6757 2.9719 2.0032 2.8600 3.6133 2.5907 2.9396 3.5362

1.7126 2.6145 3.9523 2.4856 1.5390 2.7031 2.4992 2.8452 3.6777

2.2090 3.4128 2.9869 2.3059 2.2293 3.7875 1.9334 2.9924 3.8457

5.7100 6.4276 9.7896 6.5771 6.1264 9.9584 7.0177 9.6699 8.8186

and P ij , roulette wheel selection operator presents desired number of individuals. A similar procedure as the one proposed for the NSGA-II algorithm is employed in this research to stop the NRGA algorithm. In addition, like in the NSGA-II algorithm, infeasible solutions are penalized using the death penalty approach. Fig. 4 shows the chromosomes of the NSGA-II and the NRGA algorithms, where they are presented with real numbers. In addition,

J. Sadeghi et al. / Computers & Operations Research 41 (2014) 53–64

59

Fig. 5. The mean S/N plot for different levels of the NSGA-II parameters.

Fig. 6. The mean S/N plot for different levels of the NRGA parameters.

Table 7 Tuned parameters for NSGA-II and NRGA.

5. Parameter tuning

NSGA-II

Value

NRGA

Value

Np It Pc Pm

65 1000 0.6 0.2

Np It Pc Pm

65 1000 0.6 0.2

in this research a partial-mapped crossover (PMX) operator that acts like a simple two-point crossover with a repairing procedure [36] is employed in both algorithms. Moreover, the mutation operators used in this paper are swap, reversion, and insertion, where they are randomly selected.

In order to obtain solutions with better quality (fitness function value), the parameters of both algorithms are calibrated in this section using the Taguchi method [10,37–39]. Fisher presented design of experiments called factorial designs to investigate the impact of several factors on the mean of a response [40]. Here the response is the fitness value of a solution and the factors are the parameters of the solution algorithms. Taguchi studied fractional factorial experiments (FFEs) to reduce large number of experiments in the full factorial designs [40]. In the Taguchi method, the affecting factors are divided into two parts: noise factors N and controllable factors S; where only S can be directly controlled in the experiments. Taguchi proposed a procedure to control N for reducing the

60

J. Sadeghi et al. / Computers & Operations Research 41 (2014) 53–64

Table 8 Algorithms comparison results. Problems

NSGA-II

NRGA

NPS

GD

Sp

Best_Sol

Time

Q(A, B)

NPS

GD

Sp

Best_Sol

Time

Q(B, A)

03–06 04–07 05–15 07–05 09–07 11–05 13–07 15–10 17–13 20–10

57 59 58 55 63 56 54 51 50 50

9.298 7.336 7.284 9.509 6.883 8.609 7.754 9.164 9.492 9.804

3776.545 2701.488 1943.819 4289.505 2623.713 3103.212 3718.308 3574.556 4581.100 4434.908

0.666 0.738 0.838 0.656 0.723 0.702 0.768 0.807 0.801 0.804

813.616 837.029 848.855 821.195 848.812 810.191 808.586 860.144 819.954 837.055

1 0.789 1 1 0 0 1 0 1 0

54 63 61 57 60 54 51 52 45 45

9.085 8.024 7.573 8.273 6.488 8.874 8.500 9.786 9.224 8.498

3986.192 3390.617 2779.977 3601.483 1915.397 3823.371 3740.477 3448.279 3680.501 3157.313

0.668 0.673 0.818 0.698 0.763 0.701 0.746 0.783 0.842 0.829

822.338 816.532 840.053 821.129 800.502 802.948 818.166 825.959 817.586 813.015

0 0.210 0 0 1 1 0 1 0 1

Ave.

55.3

8.513

3474.715

0.750

830.544

0.556

54.2

8.432

3352.361

0.752

817.823

0.444

variation around the target with regard to orthogonal arrays. The design that is impressed less by noise factors is the robust design. There are two ways to perform analysis of results in the Taguchi method; first is using analysis of variance for experiments with a single replicate; the second is to use the signal to noise ratio (S/N) for experiments with multiple replications. Since the one with multiple replications has a better performance, the S/N is applied in this research to analyse the solutions. For more information regarding the Taguchi method see [41].

5.1. Taguchi method implementation The parameters of both the NSGA-II and the NRGA algorithms along with their levels are shown in Table 1. The Taguchi L9 orthogonal array is utilized to run the experiments and tune the parameters. Four different responses, each representing a specific quality of a solution obtained using a MOEA, are considered for the experiments. They are (1) number of Pareto solutions (NPS), general distance of a solution (GD), spacing (Sp), and the best solution obtained (Best_Sol) [34]. Note that in order to obtain the Best_Sol, the fitness values of the solution for the two objective functions are first normalized using the linear dimensionless approach of multi-attribute decision-making processes. Then, they are summed together with a weight of 0.5 each. The solution with the highest combined value of the normalized objective function values is selected the best. Note that since a minimization problem is involved, on the one hand GD and Sp with lower values explain better efficiency of MOEAs. On the other hand, NPS and Best_Sol with bigger values explain the better efficiency of MOEAs. Regarding the basic data in a numerical example with four retailers and seven machines shown in Table 2 and each parameter level of the L9 from orthogonal arrays, each problem is run five times. Table 3 shows the results obtained based on an arbitrary run. The sums of the obtaind metrics of the five runs using the NSGA-II algorithm are normilized in Table 4. These normalized values are used in the Taguchi method as the response values of experiments based on different combinations of the parameter levels. Since a solution with the highest Sum is desired, the aim is to find maximum S/N calculated by ! 1 n 1 S=N ¼ 10 log ð30Þ ∑ n i ¼ 1 Sum2i where Sumi ; i ¼ 1; …; 5 is the response in ith replication of the Taguchi method and n ¼ 5 is the number of replications in experiments.

Tables 5 and 6 show the experimental results of NSGA-II and NRGA under different scenarios of the parameter combinations, respectively where “1”, “2”, and “3” refer to the first, the second, and the third level of the parameters. Regarding Eq. (30), these tables present S/Ns as well. In addition, it can be seen from Figs. 5 and 6 that the highest mean of S/N is the best. Therefore, Table 7 contains the optimal parameter values of the algorithms. Note that depending on the problem size and limited CPU time, Np and It may vary. Note also that Minitab 15.1.30.0. is used to employ the Taguchi method.

6. Performance comparison In order to verify and validate the results obtained by NSGA-II, its performance is compared with the one of the NRGA algorithm in this section. To do this, each problem is run 10 times, where the best run is selected for comparison. The performance measure is obtained by the simple additive weighting algorithm in multicriteria decision-making environment with equal weight of criteria defined in the previous section. In addition, the set coverage, Q(A, B), and the CPU time, Time, metrics are used to evaluate the performances as well. The set coverage metric is used to compare two Pareto optimal solutions (A, B) based on the concept of dominance relations [34]. Thus, six performance measures are used to compare NSGA-II and NRGA. Table 8 shows the results of these metrics obtained by the two algorithms using 10 different problems, where the first column represents the problem size so that for instance “03–06” refers to a problem with 3 retailers and 6 machines. The results in Table 8 show that NSGA-II performs better that NRGA in terms of the majority of metrics. The implementation results of both algorithms on the “07–05” problem are shown in Table 9, where the obtained values of both objectives are given. Fig. 7 shows these results graphically. Based on the results in Table 9 (or Fig. 7), majority of Pareto solutions obtained by NRGA are dominated by NSGA-II. 6.1. Evaluating Pareto optimal solutions of NSGA-II In order to evaluate and assess the quality of the solutions provided by NSGA-II, a genetic algorithm (GA) is employed in this section to provide lower and upper bounds resulted by solving two single-objective optimization problems, separately. The results that are shown in Table 10 show that the Pareto optimal solutions have enough adequacy such that there is only 3% difference on the average. Note that the operators used in the GA are similar to the ones used in NSGA-II of this research [8,10]. Note also that throughout this paper all calculations are performed on a PC with

J. Sadeghi et al. / Computers & Operations Research 41 (2014) 53–64

Table 9 Pareto solutions of the “07–05” problem. Solution number

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57

NSGA-II

NRGA

TC

R

TC

R

651632 203257 606460 633360 638774 489003 438417 598680 320837 433339 393122 315759 456689 411517 515903 353028 275665 353151 561288 574146 366009 483589 588188 521194 534052 406103 275542 510702 501861 424375 620502 393245 266701 306795 235448 478511 293937 371423 288523 235571 339109 470608 253843 248429 574269 553508 541955 221529 534175 555997 502922 307979 425559 398323 379203 – –

0.3111 0.0100 0.2743 0.2828 0.3037 0.1660 0.1288 0.2545 0.0547 0.1139 0.0839 0.0408 0.1417 0.1050 0.1817 0.0629 0.0305 0.0715 0.2278 0.2341 0.0737 0.1551 0.2494 0.1909 0.1962 0.0985 0.0203 0.1767 0.1706 0.1083 0.2761 0.0955 0.0196 0.0348 0.0115 0.1537 0.0336 0.0787 0.0316 0.0173 0.0602 0.1481 0.0190 0.0179 0.2385 0.2114 0.2035 0.0110 0.1999 0.2169 0.1736 0.0400 0.1119 0.0965 0.0803 – –

688787 272726 643615 430277 670515 683709 475572 638537 357992 470494 352914 430400 448672 493844 443258 558349 390183 312820 598443 408578 625343 403164 390306 507763 580171 593029 553058 665437 540077 533938 312697 461530 520744 376264 526035 290998 280506 566129 611424 611301 579110 515666 571330 345134 343950 547857 462714 526158 325678 507886 593152 303856 298778 571207 547980 336170 603521

0.3111 0.0173 0.2743 0.0839 0.2828 0.3085 0.1288 0.2562 0.0547 0.1139 0.0408 0.0955 0.1050 0.1417 0.0985 0.1909 0.0629 0.0305 0.2278 0.0787 0.2494 0.0737 0.0715 0.1481 0.2071 0.2128 0.1817 0.2805 0.1736 0.1690 0.0203 0.1083 0.1551 0.0602 0.1629 0.0190 0.0177 0.1944 0.2385 0.2341 0.2035 0.1537 0.1999 0.0400 0.0348 0.1767 0.1119 0.1660 0.0316 0.1509 0.2169 0.0196 0.0194 0.1962 0.1801 0.0340 0.2299

2.2 GHz Intel Core 2 Duo CPU, and 4 GB of RAM memory. Moreover, the algorithms were coded using the MATLAB 2010b software. 6.2. Improvement of the algorithms: hybrid NSGA-II and hybrid NRGA A popular improvement of meta-heuristics is hybridization, which mixes local search or exact algorithms and meta-heuristics for more effective and efficient problem-solving [42]. Combining a local searcher with an evolutionary algorithms, sometimes called memetic algorithms, has shown to be successful (see [43] for more

61

details). The main benefit of a local search algorithm is to prevent being trapped at a local optimal that is worse than the global one. Population-based and single-solution based meta-heuristics are the main categories in mathematical optimization to find a good solution to a problem. Since single-solution-based meta-heuristics, which sometimes called trajectory methods, begin with a single initial solution and improve it iteratively, they can be used as suitable extensions of local search algorithms such as simulated annealing (SA) [42]. This subsection uses an SA to make the proposed algorithms faster in finding near-optimum solutions. The proposed hybrid algorithms use SA to create a feasible solution.

6.2.1. SA implementation SA has its origin initiated by statistical mechanics named Metropolis algorithm [44]. Kirkpatrick et al. [45] and Černý [46] were the first who proposed SA independently. SA inspiration comes from an annealing process in metallurgy, a technique consisting of heating and controlled cooling of a material to obtain a strong crystalline structure. After carrying a material at high temperature, careful and slow cooling creates strong crystals. SA simulates these energy changes to optimize problems so that it imitates a cooling process in a system subjected to reach an equilibrium state (steady frozen state or near-optimal solution). To start the local searcher, after setting the temperature parameter ðT 0 ¼ 1000Þ; SA creates an initial solution randomly, and then iterates the following cycle until satisfying a termination criterion, here the number of iteration ðIt ¼ 10 000Þ: The operators, namely swap, reversion, and insertion, first generate a random neighbor ðsol′Þ from an initial solution ðsolÞ. Then, their cost functions denoted by f ðsol′Þ and f ðsolÞ; are obtained, respectively. If f ðsol′Þ r f ðsolÞ; sol′ is accepted as a new solution, and is replaced with sol: Otherwise, if f ðsol′Þ Zf ðsolÞ; sol′ can be accepted with the probability pobtained using the Boltzmann distribution shown in Eq. (31) [45]. Tis gradually decreased during the search process; thus, the probability of accepting inefficient solutions is high at the beginning of the search, and it gradually decreases: p ¼ expððf ðsol′Þf ðsolÞÞ=TÞ

ð31Þ

Note that the parameters of the SA are set by a trial and error method in this paper.

6.2.2. Improved results Considering SA as a local searcher, the proposed algorithms, NSGA-II and NRGA, explained in Section 4, are improved to optimize a bi-objective model presented in Eq. (27). These algorithms are called hybrid NSGA-II and hybrid NRGA thereafter. Since SA is a single-solution-based meta-heuristic, as a local searcher it can obtain feasible solution quickly. The results presented in Table 11 illustrate that SA makes the hybrid algorithms faster in terms of the required CPU time. In this case, note that the solution quality remains almost unchanged. Although the hybrid NSGA-II and hybrid NRGA become faster using the local searcher, since they have a similar cores to the ones of non-hybrid algorithms, the former algorithm performs better than the latter. Fig. 8 confirms this claim that hybrid NSGA-II is the faster algorithm. Besides, according to the results in Table 11, it can be concluded that the local searcher can reduce CPU by 0.38%. Moreover, Fig. 9 displays the considerable improvements in speed for the hybrid NSGA-II versus non-hybrid NSGA-II. Note that tuned setting of the hybrid algorithms is Np ¼50; It ¼1000; Pc ¼ 0.6; Pm ¼0.2 based on the Taguchi method.

62

J. Sadeghi et al. / Computers & Operations Research 41 (2014) 53–64

Fig. 7. Pareto solutions of the “07–05” problem.

Table 10 Evaluating Pareto optimal solutions. Problems

Best solution for the first objective

Best solution for the second objective

Optimizing by NSGA-II

Optimizing by GA

Difference rate %

Optimizing by NSGA-II

Optimizing by GA

Difference rate %

03–06 04–07 05–15 07–05 09–07 11–05 13–07 15–10 17–13 20–10

223699.83 219468.24 445098.61 203257.18 306040.16 163980.18 281314.92 474492.72 581422.96 489875.12

222773.80 221290.38 434455.84 195740.08 299015.57 161857.18 281425.92 431799.84 558422.12 465845.10

0.42 0.82 2.45 3.84 2.35 1.31 0.04 9.89 4.12 5.16

0.3113 0.2568 0.0092 0.2956 0.0448 0.5578 0.2552 0.0019 0.0223 0.0778

0.3110 0.2568 0.0096 0.3110 0.0448 0.5766 0.2576 0.0020 0.0227 0.0846

0.1 0 4.32 4.96 0 3.25 0.93 7.66 1.83 7.96

Average

338864.99

282934.91

3.04

0.1833

0.1877

3.101

Table 11 Improved results. Problems

CPU time

Best_Sol

Hybrid NSGA-II

NSGA-II

Hybrid NRGA

NRGA

Hybrid NSGA-II

NSGA-II

Hybrid NRGA

NRGA

03–06 04–07 05–15 07–05 09–07 11–05 13–07 15–10 17–13 20–10

497.00 499.30 505.09 507.62 508.45 509.15 510.64 511.82 513.49 515.07

813.62 837.03 848.86 821.20 848.81 810.19 808.59 860.14 819.95 837.06

501.00 506.30 513.09 511.62 515.45 518.15 515.64 515.82 521.49 522.07

822.34 816.53 840.05 821.13 800.50 802.95 818.17 825.96 817.59 813.02

0.666 0.739 0.836 0.654 0.723 0.703 0.766 0.809 0.822 0.829

0.666 0.738 0.838 0.656 0.723 0.702 0.768 0.807 0.801 0.804

0.668 0.673 0.818 0.697 0.764 0.701 0.745 0.784 0.842 0.83

0.668 0.673 0.818 0.698 0.763 0.701 0.746 0.783 0.842 0.829

Average

507.76

830.54

514.06

817.82

0.755

0.750

0.752

0.752

7. Discussion and conclusion In this paper, a multi-objective combinatorial optimization model of a supply chain problem including one-vendor multi-retailers considering a VMI approach was presented. The VMI model of Zavanella and Zanoni [4] was extended to include both the transportation cost and the redundancy allocation of production machines.

In the developed model, the vendor manages the retailers' inventory level and delivers the goods to the retailers in order to reduce the total inventory cost under limitations of replenishment frequency, the warehouse space for all retailers, and the average inventory level of the vendor. In order to reduce the transportation cost, the shortest route was determined in this research to deliver goods to the retailers. Inventory cost minimization was defined the

J. Sadeghi et al. / Computers & Operations Research 41 (2014) 53–64

63

Acknowledgments The authors would like to thank anonymous reviewers for their comments and suggestions that improved the presentation of the paper. References

Fig. 8. Boxplot of hybrid NSGA-II and hybrid NRGA.

Fig. 9. Boxplot of hybrid NSGA-II and NSGA-II.

first objective, and maximization of the system reliability of the machines that produce the good was defined the second. Since the developed model of the problem obtained in (27) was NP-hard, two multi-objective genetic algorithms, namely NSGA-II and NRGA, were used to find Pareto fronts. After calibrating the parameters of both algorithms using the Taguchi method, the algorithms were compared in terms of six performance metrics. Based on the results in Table 8, the two algorithms had similar performances with 3.5% tolerance in solutions. However, NSGA-II had 25% reduction in the mean of the set coverage metric. In other words, NSGA-II dominates majority of the Pareto optimal solutions of NRGA. Moreover, a GA was employed to assess Pareto optimal solutions of NSGAII. Besides, a local searcher, SA, was integrated to both algorithms to make them hybrid. The results of employing the hybrid algorithms showed while the CPU time to find a near-optimum solution would become significantly less, the performances in terms of the “bestsolution” obtained remain almost unchanged. Some future research in this area are recommend as follows: i. Developing a model for non-deterministic condition, such as stochastic or fuzzy demand. ii. Considering the VMI model for several warehouses and multiple vendors. iii. Applying response surface methodology (RSM) to tune the parameters. iv. Considering shortages and discounts in the modeling.

[1] Hugos M. Essentials of supply chain management. 2 ed. New Jersey: Wiley; 2006. [2] Simchi-Levi D, Kaminsky P, Simchi-Levi E, Shankar R. Designing and managing the supply chain: concepts, strategies and case studies. 3rd ed. New York: McGraw-Hill; 2007. [3] Disney SM, Towill DR. A procedure for the optimization of the dynamic response of a vendor managed inventory system. Computers & Industrial Engineering 2002;43(1–2):27–58. [4] Zavanella L, Zanoni S. A one-vendor multi-buyer integrated productioninventory model: the consignment stock case. International Journal of Production Economics 2009;118(1):225–32. [5] Kuo W, Zuo MJ. Optimal reliability modeling: principles and applications.New Jersey: Wiley; 2002. [6] Yao Y, Evers PT, Dresner ME. Supply chain integration in vendor-managed inventory. Decision Support Systems 2007;43(2):663–74. [7] Darwish MA, Odah OM. Vendor managed inventory model for single-vendor multi-retailer supply chains. European Journal of Operational Research 2010;204(3):473–84. [8] Pasandideh SHR, Niaki STA, Nia AR. A genetic algorithm for vendor managed inventory control system of multi-product multi-constraint economic order quantity model. Expert Systems with Applications 2011;38(3):2708–16. [9] Cárdenas-Barrón LE, Treviño-Garza G, Wee HM. A simple and better algorithm to solve the vendor managed inventory control system of multi-product multiconstraint economic order quantity model. Expert Systems with Applications 2012;39(3):3888–95. [10] Sadeghi J, Sadeghi A, Saidi-Mehrabad M. A parameter-tuned genetic algorithm for vendor managed inventory model for a case single-vendor single-retailer with multi-product and multi-constraint. Journal of Optimization in Industrial Engineering 2011;4(9):57–67. [11] Hariga M, Gumus M, Daghfous A, Goyal SK. A vendor managed inventory model under contractual storage agreement. Computers & Operations Research 2013;40(8):2138–44. [12] Kiesmüller GP. Broekmeulen RACM. The benefit of VMI strategies in a stochastic multi-product serial two echelon system. Computers & Operations Research 2010;37(2):406–16. [13] Goyal SK. A joint economic-lot-size model for purchaser and vendor: a comment. Decision Sciences 1988;19(1):236–41. [14] Goyal SK. A one-vendor multi-buyer integrated inventory model: A comment. European Journal of Operational Research 1995;82(1):209–10. [15] Lu L. A one-vendor multi-buyer integrated inventory model. European Journal of Operational Research 1995;81(2):312–23. [16] Hill RM. The optimal production and shipment policy for the single-vendor single-buyer integrated production-inventory model. International Journal of Production Research 1999;37(11):2463–75. [17] Hill RM, Omar M. Another look at the single-vendor single-buyer integrated production-inventory problem. International Journal of Production Research 2006;44(4):791–800. [18] Braglia M, Zavanella L. Modelling an industrial strategy for inventory management in supply chains: the ‘consignment stock' case. International Journal of Production Research 2003;41(16):3793–808. [19] Nachiappan SP, Jawahar N. A genetic algorithm for optimal operating parameters of VMI system in a two-echelon supply chain. European Journal of Operational Research 2007;182(3):1433–52. [20] Michaelraj LA, Shahabudeen P. Replenishment policies for sustainable business development in a continuous credit based vendor managed inventory distribution system. Computers & Industrial Engineering 2009;56(1):260–6. [21] Gupta RK, Bhunia AK, Goyal SK. An application of genetic algorithm in solving an inventory model with advance payment and interval valued inventory costs. Mathematical and Computer Modelling 2009;49(5–6):893–905. [22] Yang W, Chan FTS, Kumar V. Optimizing replenishment polices using genetic algorithm for single-warehouse multi-retailer system. Expert Systems with Applications 2012;39(3):3081–6. [23] Yu Y, Huang GQ, Hong Z, Zhang X. An integrated pricing and deteriorating model and a hybrid algorithm for a VMI (vendor-managed-inventory) supply chain. IEEE Transactions on Automation Science and Engineering 2011;8 (4):673–82. [24] Yu Y, Hong Z, Zhang LL, Liang L, Chu C. Optimal selection of retailers for a manufacturing vendor in a vendor managed inventory system. European Journal of Operational Research 2013;225(2):273–84. [25] Kuo RJ, Han YS. A hybrid of genetic algorithm and particle swarm optimization for solving bi-level linear programming problem – a case study on supply chain model. Applied Mathematical Modelling 2011;35(8):3905–17. [26] Liao S-H, Hsieh C-L, Lai P-J. An evolutionary approach for multi-objective optimization of the integrated location-inventory distribution network

64

[27] [28] [29] [30] [31] [32]

[33]

[34] [35]

[36]

J. Sadeghi et al. / Computers & Operations Research 41 (2014) 53–64

problem in vendor-managed inventory. Expert Systems with Applications 2011;38(6):6768–76. Pham HE. Handbook of reliability engineering.London: Springer-Verlag London Ltd; 2003. Talbi E-G. Metaheuristics: from design to implementation. New York: Wiley; 2009. Miller RRE, Thatcher JW. Complexity of computer computations: proceedings. USA: Springer; 1972. Lin M-S, Chen D-J. The computational complexity of the reliability problem on distributed systems. Information Processing Letters 1997;64(3):143–7. Chern M-S. On the computational complexity of reliability redundancy allocation in a series system. Operations Research Letters 1992;11(5):309–15. Peng T, Jian M, Dong-Mo Z. Non-linear integer programming by Darwin and Boltzmann mixed strategy. European Journal of Operational Research 1998;105(1):224–35. Deb K, Agrawal S, Pratap A, Meyarivan T. A fast elitist non-dominated sorting genetic algorithm for multi-objective optimisation: NSGA-II. In: Proceedings of the 6th international conference on parallel problem solving from nature. Springer-Verlag; 2000. p. 849–58. Deb K. Multi-objective optimization using evolutionary algorithms. 1st ed. New York: Wiley; 2001. Jadaan OA, Rajamani L, Rao CR. Non-dominated ranked genetic algorithm for solving multi-objective optimization problems: NRGA. Journal of Theoretical and Applied Information Technology 2008;2:60–7. Gen M, Cheng R. Genetic algorithms and engineering design. 1st ed.New York: Wiley; 2000.

[37] Najafi AA, Niaki STA, Shahsavar M. A parameter-tuned genetic algorithm for the resource investment problem with discounted cash flows and generalized precedence relations. Computers & Operations Research 2009;36 (11):2994–3001. [38] Sels V, Vanhoucke M. A hybrid genetic algorithm for the single machine maximum lateness problem with release times and family setups. Computers & Operations Research 2012;39(10):2346–58. [39] Fazel Zarandi MH, Mosadegh H, Fattahi M. Two-machine robotic cell scheduling problem with sequence-dependent setup times. Computers & Operations Research 2013;40(5):1420–34. [40] Roy R. A primer on the Taguchi method. New York: Society of Manufacturing Engineers; 1990. [41] Taguchi G, Chowdhury S, Wu Y. Taguchi's quality engineering handbook.New Jersey: Wiley; 2005. [42] Boussaïd I, Lepagnot J, Siarry P. A survey on optimization metaheuristics. Information Sciences 2013;237(0):82–117. [43] Blum C, Puchinger J, Raidl GR, Roli A. Hybrid metaheuristics in combinatorial optimization: a survey. Applied Soft Computing 2011;11(6):4135–51. [44] Metropolis N, Rosenbluth AW, Rosenbluth MN, Teller AH, Teller E. Equation of state calculations by fast computing machines. The Journal of Chemical Physics 1953;21(6):1087–92. [45] Kirkpatrick S, Gelatt CD, Vecchi MP. Optimization by simulated annealing. Science 1983;220(4598):671–80. [46] Černý V. Thermodynamical approach to the traveling salesman problem: an efficient simulation algorithm. Journal of Optimization Theory and Applications 1985;45(1):41–51.