Modeling and Optimization of Stochastic Joint Replenishment and Delivery Scheduling Problem with Uncertain Costs

The stochastic joint replenishment and delivery scheduling (JRD) problem is a key issue in supply chainmanagement and is amajor concern for companies. So far, all of the work on stochastic JRDs is under explicit environment. However, the decisionmakers often have to face vague operational conditions.We develop a practical JRDmodel with stochastic demand under fuzzy backlogging cost, fuzzy minor ordering cost, and fuzzy inventory holding cost. The problem is to determine procedures for inventory management and vehicle routing simultaneously so that the warehouse may satisfy demand at a minimum long-run average cost. Subsequently, the fuzzy total cost is defuzzified by the gradedmean integration representation and centroid approaches to rank fuzzy numbers. To find optimal coordinated decisions, a modified adaptive differential evolution algorithm (MADE) is utilized to find the minimum long-run average total cost. Results of numerical examples indicate that the proposed JRD model can be used to simulate fuzzy environment efficiently, and the MADE outperforms genetic algorithm with a lower total cost and higher convergence rate. The proposed methods can be applied to many industries and can help obtaining optimal decisions under uncertain environment.


Introduction
The joint replenishment problem (JRP) has been heavily researched since the early work of Shu [1].The JRP means to group items into the same order from a supplier to achieve the purpose of sharing the main preparation costs and saving the procurement costs (Goyal [2]).The existing research of JRPs can be classified into two categories: (1) the classic JRPs under constant demand, (2) JRPs under stochastic or dynamic demand.Corresponding review is available in Khouja and Goyal [3] and Robinson et al. [4].
Supply chain management has received much attention and global sourcing has been widespread in recent years.Since cost savings can be achieved by using joint replenishment policy, managers have realized that jointly considering JRP and delivery scheduling (JRD) can obtain a scale effect of replenishment and transportation simultaneously, thus further reducing total cost (Qu et al. [5], Wang et al. [6]).Unfortunately, the literature of JRD is limited and most of scholars studied JRDs with deterministic demand.Sindhuchao et al. [7] developed the inventory and transportation policy for the central warehouse and decentralized suppliers.They considered the transportation capacity and transportation frequency constraints.Chan et al. [8] discussed delivery scheduling after the JRP has been solved and, respectively, designed four different objectives.Cha et al. [9] developed a joint replenishment and transportation model in a onewarehouse, n-retailer system based on improved well-known heuristic named RAND and Genetic Algorithm (GA).They also extended their model with constraints and showed the flexibility of GA.Moon et al. [10] proposed two algorithms to determine the joint replenishment and delivery policies from the perspective of a third party warehouse.Unlike the studies above that consider deterministic demand, Qu et al. [5] discussed a multiitem JRD with modified periodic-review policy under stochastic demand.They considered nonlinear transportation cost and proposed an efficient heuristic algorithm to solve the problem.However, all the existing JRD models mentioned above assume explicit vital factors, which is too harsh in practice.
In reality, managers often have to make decisions under imprecise operational conditions (Zeng et al. [11]).Uncertainties have been restrained in a variety of ways and fuzzy set theory has a long history for managing inventories.Chen and Chang [12] proposed a fuzzy economic production quantity model with defective productions that cannot be repaired.Wang et al. [13] studied continuous review inventory models in which a fraction of demand is backordered and the remaining fraction is lost during the stock-out period under fuzzy demands.Dey and Chakraborty [14] discussed a single period inventory model under the mixed fuzzy random environment and provided a methodology to handle the model.On the other hand, defuzzification has been a favorite approach in many inventory studies for its simplicity.Defuzzification can easily transfer fuzziness to be explicit without complex analysis (Lin [15]).So, we will also utilize two defuzzification approaches to handle the fuzziness.However, to the best of our knowledge, only a few papers discussed JRD under fuzzy environment (Wang et al. [16]).However, they assume that demand rates are known and constant.It is too ideal and results in the relatively low practical value of the proposed model.
Moreover, it is very difficult to find optimal solutions for the defuzzified model effectively.The JRP has been proven to be an NP-hard problem (Arkin et al., [17]).While optimal solutions are very hard to be obtained, many acceptable solutions are designed in practice.Current approaches for solving JRPs include an iterative algorithm (Goyal [18]), power-oftwo (PoT) policy (Lee and Yao [19]), RAND algorithm (Kaspi and Rosenblatt [20]), evolutionary algorithm (Olsen [21]), branch-and-bound (Robinson and Gao [22]), and simulated annealing metaheuristic (Robinson et al. [23]).For the JRD model, there are more decision variables which will influence each other.Unfortunately, defuzzification also increases the complexity of the problem further.It is more difficult to solve the fuzzy JRD model in generally acceptable time and CPU memory requirements.Several heuristics may give acceptable solutions when the scale of optimization is small.However, it is not easy to find a proper heuristic with a robust performance.
On the other hand, metaheuristics have grown quickly in order to attain better solutions.Among these algorithms, the evolutionary algorithms (EAs) and especially genetic algorithm (GA) have been proved to be effective algorithms for the JRD.For example, Cha et al. [9] designed two heuristics and a hybrid GA.Results showed the GA outperformed the best solution of each problem in their experiment (the maximum error found is 0.3660% and the average error found is only 0.0268%).However, the GA displays inherent difficulties in performing local search for some numerical applications and the optimal solution cannot be obtained.So, it is necessary to find a novel algorithm to deal with the fuzzy JRD more effectively.
The aim of this paper is to model and optimize the practical JRD policy with stochastic demands under uncertainty.This topic is interesting because of the widely adoption of the JRD policy in many industries and the operability and rationality of the method to handle uncertain costs.This study makes the work of Qu et al. [5] and Wang et al. [6] become more practical by considering the fuzziness of minor ordering cost, holding cost, and backlogging cost.This study also extends the work of Wang et al. [16] which studied fuzzy JRD under determined environment.The introduction of stochastic demands in the JRD results in rather complex mathematical properties of the proposal model.So, a new approach should be designed to solve this NP-hard problem better.Therefore, we propose a simple and effective modified adaptive differential evolution algorithm (MADE) for the fuzzy JRD model firstly.Results of an example show that an MADE performs better than GA which is also regarded as an effective algorithm for this problem.Then, the MADE is utilized to solve the fuzzy JRD and the results of fuzzy JRD are compared with the JRD.Numerical examples show that the ranges of fuzzy parameters will influence the decision of this problem.
The rest of this paper is organized as follows.In Section 2, we propose a fuzzy JRD model and give the results by two defuzzification approaches.Section 3 proposes an MADE to solve the defuzzified JRD model.Section 4 presents numerical studies and analysis.Section 5 contains the conclusion and provides directions for future research.Preliminaries are presented in the appendix.  : the order interval of the item  (decision variable).

Fuzzy JRD Model and Analysis
Sometimes, the stock-out factor   is predetermined by decision maker.However, in order to obtain the optimal value of   or service level, Eynan and Kropp [25] regarded   as a decision variable, which is the same with the study of Qu et al. [5].In this study, we also regard   as a decision variable.

Formulation of JRD.
The total cost consists of the following parts: (1) the inventory cost, which includes ordering, holding, and backlog cost; (2) the transportation cost, which includes dispatching and stopover and routing costs.The details of different costs will be discussed as follows.
(1) Inventory Cost.The ordering cost consists of the major ordering cost which is independent of the number of different items in the order and the minor ordering cost which depends on the number of different items in the order.Hence, the annual ordering cost is We suppose the demand of each item follows the normal distribution and the lead time is   .As described by Eynan and Kropp [25] and Qu et al. [5], when decision variables   and  are determined, the demand of item  will be normally distributed in each replenishment interval (   +   ), and the expectation  is equal to   (   +   ); the variance Var is equal to   (   +   ).Take into the lead time, the maximum inventory   in each replenishment interval of item  includes the expected demand during the time and the safety stock expressed as a multiple of the standard deviation of demand forecast errors during period (   +   ), that is to say   =   (   +   ) +   √  (   +   ).The average net inventory level of item  during the interval is approximately   −   (  +   /2) = (1/2)     +   √  ( +   ); thus, the average inventory holding cost for periodic policy is A stockout occurs following an order placed at time  if the cumulative demand   (between  and  +   +   ) exceeds   .Then, there will be a penalty cost   for each unit backlogged.If demand has a density function  (  ,   +   ) over the interval of length   +  , the mean backlogging cost in the long run is (2) Transportation Cost.Define  as the lowest common multiple of   .All items will repeat the same replenishment and distribution policy in every other regeneration cycle .So we only need to consider the replenishment and distribution policy in a regeneration cycle .
For each replenishment period , the transportation cost includes the costs of stopover (at those plants visited) plus the travelling costs.Over a regeneration cycle, the average transportation cost is where 1, if the items are kept in the supplier  in the th replenishment period 0, otherwise. ( (3) The Objective.Then, the objective is to minimize the total long-run expected cost (TC); where

Fuzzy JRD Involving Fuzzy Minor Ordering Cost, Holding
Cost, and Backlogging Cost 2.2.1.Formulation.As reported in Vujošević et al. [26], an interesting real problem is when human originated data like holding cost and ordering cost are not precisely known but subjectively estimated or linguistically expressed because of the lack of the accurate history data.Suppose that the minor ordering cost   and holding cost ℎ  are not precisely known.For example, let the linguistic estimates of   and ℎ  be as follows.
"The minor ordering cost is about 8 dollars per order, but not less than 6 dollars per order nor more than 15 dollars per order." "The holding cost is between 15% and 20% of the average inventory value, but not less than 12% nor more than 25%." These imprecise definements in linguistics are always described by fuzzy numbers.So, fuzzy variables are also utilized to handle the JRD problem under uncertainty.Two popular kinds of fuzzy numbers characterized by triangular and trapezoidal membership functions (MFs) are presented in the appendix.They were widely used to solve uncertainty problem because of their intuitive appeal and their perceived computational efficacy (Pramanik and Biswas [27]; Giachetti and Young [28]).In addition, trapezoidal fuzzy numbers can be obtained by extending triangular fuzzy numbers (TFNs).That is to say, TFNs can be considered as a particular kind of trapezoidal fuzzy numbers (Vujošević et al. [26]; Giachetti and Young [28]).Thus, we treat the minor ordering cost   , holding cost ℎ  , and backlogging cost   as trapezoidal fuzzy numbers in this study; that is, According to Definition A.4 in the appendix, the MFs of s , h , and τ are as follows: Take each TFN into (6) to get the fuzzy number of TC.For the simplify, we denote that Hence, 2.2.2.Defuzzified TC.Defuzzification has been a favorite approach in many inventory studies for its simplicity.Defuzzification can easily transfer fuzziness to be explicit without complex analysis.In this study, two approaches will be utilized to defuzzify the fuzzy total cost.
(1) Graded Mean Integration Representation (GMIR) Approach.Chen and Hsieh [29] introduced the GMIR method based on the integral value of the graded mean h-level of the generalized fuzzy number for defuzzifying a generalized fuzzy number.This method is reasonable because it adopts a grade as the important degree of each point of the support set of a fuzzy number for representing the fuzzy number (G. C. Mahata and P. Mahata [30]).So, we also use the GMIR approach to defuzzify the fuzzy number.
where ℎ is an arbitrary value in the interval of the membership function.According to ( 12) and ( 11), we can obtain the defuzzified total cost So, the simplified objective is (2) Centroid.Centroid is a widely used approach for ranking fuzzy numbers by distance method based on calculating the centroid point, where distance means from original point to the centroid point (Wang et al. [31]).
then the defuzzified function of Ã is Hence, the simplified objective is expressed as ] . (16)

MADE Algorithm for Fuzzy JRD Model
3.1.Rationales of Using Differential Evolution.Differential evolution algorithm (DE) was proposed by Storn and Price [32] for complex continuous nonlinear, nondifferentiable, and multimodal optimization problem.This technique combines simple arithmetic operators with the classical events of crossover, mutation, and selection to evolve from a randomly generated starting population to a final solution.Due to its simple structure, easy implementation, quick convergence, and robustness, DE has been turned out to be one of the best evolution algorithms in a variety of fields (Wang et al. [33]; Wee et al. [34]; Zheng and Yamashiro [35]).However, the mutation factor  in a typical DE is a constant, which is actually difficult to determine.The main role of  is to perturb the individual in order to avoid being trapped in a local optimum.So, if  is too small, the degree of population diversity cannot reach a certain value and it is easy to be premature.Vice versa, if  is too large, the chance to find optimal solution will decrease, and the effectiveness of the algorithm will be very low in large decision space.More description about the value of  effects on evolution algorithm can be found in Storn and Price [32] and Salman et al. [36].In this study, we use a modified adaptive DE to improve the quality of solutions of the defuzzified JRD.According to our rule, in the early stage of search,  is large enough to keep the population diversity and avoid premature.Then,  is decreased gradually, and in the latter stage  is small enough to remain the excellent individual and to increase the chance to converge to the optimal solution.

MADE-Based
Procedures for Fuzzy JRD.The model aforementioned in Section 2.2.2 can be decomposed into a JRP and a transportation problem.We use MADE to determine the optimal combination of   ,   , and , and then we can obtain a satisfied solution of the annual total cost.The detailed steps are as follows: Step 1. Representation and initialization.Let   represent the population size of a target population, and let  , denote the th target individual at iteration .
Because   is an integer, we set the lower bound  LB  = 1, which are clearly the lower bounds of item .For the upper bound  UB  , we use values of 100 to ensure that the optimal solution will not escape according to the related experience of Cha et al. [9].According to the properties of the standard normal distribution, we have ( ≤   ) = 1 − ( >   ) = 0.995 when   = 2.575.So, the range of   is set as [ LB , UB ] = [0, 4] which will cover more than 99.99% of the demand.
According to the related experience of Khouja and Goyal [3] and Qu et al. [5],  is ranged from 0 to 1.It is a fact that the annual demand of each item   is larger than 600 in the following examples; so it is reasonable to assume that at least one item will be replenished per year when   is not less than 1.That is to say,  should be less than 1.
By combining   ,   , and , we can obtain the th individual: Then, create the initial population randomly.
Step 2. If the maximum number of iterations is reached, the algorithm will be stopped and output the optimal solution; otherwise, go to the next step.
Step 3 (mutation).For each target individual  , ( = 1, 2, . . .,   ), a mutant individual V ,+1 is generated according to With randomly chosen integer indexes  1 ,  2 ,  3 ∈ {1, 2, . . .,   }.Note that indexes have to be different from each other and from the running index. is called the mutation factor which controls the amplification of the differential variations ( 2 ,  ,  3 ,  ) and  can be obtained by where  min is the minimum of , and  max is the maximum of .The way to confirm these parameters will be discussed in Section 4.1.
Step 4 (crossover).The crossover operator implements a discrete recombination of the trial individual V ,+1 and the parent individual  , to produce offspring  ,+1 .The crossover is implemented as follows: where   ,  refers to the th element of the individual  , . ,+1 and V ,+1 are similarly defined.rand() is the th evaluation of a uniform random number generator between [0, 1].rnb () is a randomly chosen index from 1, 2, . . .,  which ensures that  ,+1 gets at least one parameter from V ,+1 .Otherwise, no new parent individual would be produced and the population would not alter.CR is the crossover rate between [0, 1] which must be determined by the user.
Step 5 (selection).The evaluation function of an offspring is one-to-one competition in the MADE.It means that the resulting trial individual will only replace the original if it has a lower objective function value.Otherwise, the parent will remain in the next generation.The rule is as follows: where () means ( 14) for the GMIR approach or ( 16) for the centroid approach.This is the same for all variants of the MADE.Then, the best individual of the next generation will be at least as fit as the best individual of the current generation.
Then, the procedure will return to Step 2 after setting  =  + 1.
The steps can be described by the flow chart as shown in Figure 1.be produced by (20) Mutation: the mutant individual can be generated by (18) Initialization: (k i , z i , T) Given parameters:

An Example of the Implementation of Operations of MADE.
In this section, a simple example of four-product and three-supplier problem is presented in order to better understand the implementation of MADE.
And rnb () = 6 for this target vector.Therefore, the trial vector  ,+1 is 4.2 4 1 3 3.22 1.91 1.39 1.83 0.109.Since "4.2" is not an integer, we round it to the nearest integer; that is, it will be changed to "4. " Meanwhile, "3.22" exceeds the range of variety; it is substituted by a number randomly generated from [0, 3].After that, the sequences of adjusted  ,+1 are:

Numerical Examples and Analysis
In Section 4.1, we will test the performance of MADE with results of GA which was tested to be one of the effective algorithms by a stochastic JRD.Then, the MADE will be used to handle the fuzzy JRD.

An Example of Stochastic JRD
(1) Basic Data.The data of the example is the same as Qu et al. [5].We suppose that the central warehouse is in charge of the replenishment of 4 items and the inventory management.There are 3 suppliers and the third supplier produces the third and fourth items.The major ordering cost  is 100, the unit distance cost  is 0.5, and other data are shown in Tables 1, 2, and 3, respectively.
(2) Parameters Setting and Results.Storn and Price [32] suggested that  = 0.5, CR = 0.1 is usually a good initial choice and   is appropriate for ten times the dimensionality of the problem.Some scholars also provided useful suggestions.
According to the recommendation of Neri and Tirronen [37] and Liu and Lampinen [38], we set the relevant parameters of the MADE as follows:   = 100,  min = 0.2,  max = 1.2, CR = 0.1, and Gen  = 150.Table 4 shows the result of Qu et al. [5] in comparison to our result solved by MADE and GA.We programmed the approach of [5], MADE, and GA with the language of Matlab 6.5.Comparison of the CPU time in a computer (memory: 4G; CPU: Intel Core i5-2450 M 2.50 GHz) by running the example 50 times is also shown in Table 4.
Both GA and MADE are intelligence EAs.For a perspicuous and direct understanding, we give two convergence curves as shown in Figure 2.
Results listed in Table 4 indicate that the total cost gained by MADE is the smallest and the total cost gotten by GA is the worst.We can also conclude the following: (1) the holding cost of item 1 is 5.6 which is the lowest among all items so that the solved optimal value of  is 2 and the replenishment cycle is the longest; (2) the demand, inventory cost, and shortage cost of item 3 are the highest among all items so that we should replenish it in every replenishment cycle.Obviously, MADE is better than GA and the algorithm of Qu et al. [5].So, the MADE will be used to solve the fuzzy JRD.

Three Scenarios of Fuzzy JRD Model and Analysis.
In this section, we will design three different scenarios to compare the results of JRD and fuzzy JRD.

Scenarios 1.
We set the range of the fuzzy set to where X is the value under certainty.The related data is listed in Table 5.
The major ordering cost  was set to 100, 300, 500, and 700, respectively.The results computed by the MADE under uncertainty and determinacy are shown in Table 6.Here, the average CPU time by running each example 20 times in the abovementioned computer is also reported in Table 6.Statistical test of TC is not provided because MADE is quite effective to find the best TC every time for this problem (maximum error is less than 0.002% and is regarded as an optimal solution).In the following Scenarios 2 and 3, the CPU time is obtained by the same method and the statistic of TC is also not given for the same reason.
Table 6 shows the following: (1) results of TC under determinacy are larger than TC under uncertainty; (2)   and  obtained under uncertainty are larger than the corresponding values under determinacy; (3) the best results of TC using GMIR are larger than the corresponding values using centroid; (4) results of   and  obtained by centroid are larger than GMIR; (5) the CPU times of JRD and fuzzy JRD with GMIR are very close, but the CPU times of fuzzy JRD with centroid are much larger because of the relatively complex computation process.7 shows the data of the ordering cost s , holding cost h , and shortage cost τ under uncertainty.

Scenarios 2 (Symmetrical Fuzzy MFs). Set the range of the fuzzy set
When  was set to 100, 300, 500, and 700, respectively, the results computed by the MADE under uncertainty and determinacy are shown in Table 8.
Results shown in Table 8 indicate that the results of the best TC and  under uncertainty and determinacy are the same for this situation.However, in reality, we are always not able to acquire the accurate value of the parameters or we only can acquire it with very high cost.So, we can use fuzzy theory to solve it reasonably.4.2.3.Scenarios 3. We extend the range of the fuzzy set to [d 1 X, d 2 X, d 3 X, d 4 X] = [0.9X, 0.95 X, 1.05 X, 1.3X] and the related data is shown in Table 9.
We set the major ordering cost  to 100, 300, 500, and 700, respectively, and solve the problem in the same way.Table 10 shows the comparison of results under uncertainty and determinacy after extension.
Table 10 indicates that the results of the best TC under uncertainty are larger than the corresponding results of determinacy situation.In comparison with results of  under determinacy, the basic replenishment cycle is smaller.
From Tables 6, 8, and 10, we can reasonably come to the following conclusions: (1) the ranges of fuzzy parameters will inevitably influence the decision of the JRD.So, it is important to utilize every useful information and to correct judgment to confirm the ranges of fuzzy numbers scientifically and reasonably; (2) MADE is quite effective for this fuzzy JRD problem because it can find the best minimum long-run average total cost with 100% convergence rate.The CPU time is also accepted for decision makers.

Conclusion and Future Research
This paper is an interdisciplinary research of the fuzzy inventory model and intelligent optimization algorithm.A fuzzy JRD model for the one-warehouse, n-supplier system is studied while assuming warehouse's minor ordering cost, warehouse's inventory holding cost, and warehouse's storage cost to be fuzzy numbers and a new DE is proposed.The main contributions are as follows: (1) Due to the nonavailability of sufficient and precise input data, accurate predicted values of the minor ordering cost, inventory holding cost, and backlogging cost cannot be obtained easily; while fuzzy numbers can efficiently model the imprecise values.Obviously, it is more reasonable to handle imprecise values using the fuzzy theory in the JRD model.It is the first time to introduce fuzziness into the stochastic JRD, which will widen the application field of fuzzy theory and will make the JRD become more practical.
(2) Since the stochastic JRD policy is widely used in many industries, such as manufacturing, wholesale, maintenance, repair, and operating (MRO) supplies, the proposed fuzzy JRD models with two defuzzification methods can also be applied in these industries.(3) The proposed simple and effective MADE-based approach can find the optimal cycle time and safety factor of each item for the defuzzified JRD effectively.Moreover, the convergence rate of MADE outperforms another popular GA-based approach.Therefore, managers can always know the corresponding minimum total cost under fuzzy environment, just like the case of precise data.This study expands the application field of the DE.
In the future, a dependent-chance programming model or a chance-constrained programming model can be designed for the JRD model of the one-warehouse, n-retailer system in fuzzy environment.The DE still can provide good solutions to these problems.
Selection: using (21) to choose better individual to next generation Crossover: the trail individual can

Figure 1 :
Figure 1: Flow chart of the MADE.

4 Figure 2 :
Figure 2: Convergence curve of MADE and GA.

Table 1 :
The basic data under deterministic situation.

Table 2 :
The items supplied by suppliers.

Table 3 :
The distance between the warehouse and suppliers.

Table 4 :
Results comparison of three approaches.

Table 5 :
The basic data after extension.

Table 7 :
The basic data under uncertainty.

Table 10 :
Results comparison under different situations after extension (scenarios 3).