- Email: [email protected]
The open vehicle routing problem with fuzzy demands
Expert Systems with Applications 37 (2010) 2405–2411
Contents lists available at ScienceDirect
Expert Systems with Applications journal homepage: www.elsevier.com/locate/eswa
The open vehicle routing problem with fuzzy demands Cao Erbao *, Lai Mingyong College of Economics and Trade, Hunan University, Changsha 410079, China
a r t i c l e
i n f o
Keywords: Logistics distribution Open vehicle routing problem Fuzzy credibility Stochastic simulation Differential evolution algorithm Optimization
a b s t r a c t According to the open vehicle routing problem (OVRP), a vehicle is not required to return to the distribution depot after servicing the last customer on its route. In this paper, the open vehicle routing problem with fuzzy demands (OVRPFD) is considered. A fuzzy chance-constrained program model is designed based on fuzzy credibility theory. Stochastic simulation and an improved differential evolution algorithm are integrated so as to use a hybrid intelligent algorithm to solve the OVRPFD model. The influence of the decision-maker’s preference on the final outcome of the problem is analyzed using stochastic simulation, and the range of possible preferences is calculated. Crown Copyright Ó 2009 Published by Elsevier Ltd. All rights reserved.
1. Introduction The open vehicle routing problem (OVRP) analyzes efficient routes with minimum total cost for a fleet of vehicles that serve some commodity to a given number of customers. Each customer is visited exactly once by one vehicle, while vehicle activity is bounded by capacity constraints, duration constraints and time constraints. Either each route is a sequence of customers that starts at the depot and finishes at one of the customers to whom goods are delivered, or each route is a sequence of customers that begins at a certain customer and ends at the distribution depot, where goods are gathered. The OVRP differs from the well-known vehicle routing problem (VRP) (Dantzing & Ramser, 1959) in that the vehicles do not necessarily return to their original locations after delivering goods to customers; if they do, they must visit the same customers in the reverse order. The OVRP is encountered in practice in many contexts, such as home delivery of packages and newspapers. Contractors who are not employees of the delivery company use their own vehicles and do not return to the depot; the ‘‘missing” node on a route can be driver’s home or a parking lot where the vehicle stays overnight. The major difference in theory between the OVRP and the VRP is that the routes in the OVRP consist of Hamiltonian paths originating at the depot and terminating at one of the customers, while the routes in the VRP are Hamiltonian cycles. In other words, the best Hamiltonian path is NP-hard, since the Hamiltonian path problem is equivalent to the traveling salesperson problem, which is known to be NP-hard (Reinelt, 1991). The best Hamiltonian path problem with a fixed source node must be solved for each vehicle in the OVRP, and OVRP solutions involve finding the best Hamilto* Corresponding author. Tel.: +86 7318646317; fax: +86 7318684825. E-mail address: [email protected] (E. Cao).
nian path for each set of customers assigned to a vehicle. Consequently, the OVRP is also a NP-hard problem. The OVRP has received sparse attention in the literature compared to the VRP. While the earliest description of OVRP offered by Schrage (1981) appeared in the literature over 20 years ago, OVRPs have just recently attracted the attention of practitioners and researchers. From the early 1980s to the late 1990s, the OVRP received very little attention in the operations research literature. However, since 2000, several researchers have used various heuristics and meta-heuristics to solve the OVRP with some success. Brandão (2004) and Fu, Eglese, and Li (2005) implemented a tabu search (TS) heuristic to solve the OVRP with constraints on vehicle capacity and maximum route length. Tarantilis and Kiranoudis (2002) considered a real-life fresh meat distribution problem, formulating it as a multi-depot OVRP; the problem was solved using a metaheuristic method called the list-based threshold-accepting (LBTA) algorithm. Tarantilis, Ioannou, Kiranoudis, and Prastacos (2005) solved the OVRP by adopting the LBTA algorithm they proposed for the solution of the multi-depot OVRP (Tarantilis & Kiranoudis, 2002). Pisinger and Ropke (2007) presented an adaptive largeneighborhood search heuristic to solve an OVRP. Li, Golden, and Wasil (2007) proposed a record-to-record travel heuristic and a deterministic variant of simulated annealing to solve the OVRP. However, traditional studies of the OVRP as well as the VRP have assumed that all information is deterministic, including customer information, vehicles information, roads information, and dispatcher information. As such, the proposed algorithms were only used to solve deterministic circumstances. That is, all the above-mentioned papers assumed that the demands of all customers visited on its route by any vehicle were deterministic. Actually, in some new systems, it is hard to describe the parameters of the VRP as deterministic, because there exists much
0957-4174/$ - see front matter Crown Copyright Ó 2009 Published by Elsevier Ltd. All rights reserved. doi:10.1016/j.eswa.2009.07.021
2406
E. Cao, M. Lai / Expert Systems with Applications 37 (2010) 2405–2411
uncertainty in data such as customer demands, traveling time and the set of customers to be visited. We call these problems nondeterministic VRPs. Given this aspect of VRPs, a consideration of stochastic vehicle routing problems (SVRP) and fuzzy vehicle routing problems (FVRP) may be useful. SVRPs arise whenever some elements of a given problem are random. Common examples include stochastic demands and stochastic travel times. Sometimes, the set of customers to be visited is not known with certainty. In this case, each customer has a certain probability of being visited. Researchers have developed many models and algorithms for SVRPs (Bertsimas, 1992; Dror, Laporte, & Trudeau, 1989; Gendreau, Laporate, & Seguin, 1996; Liu & Lai, 2002). Alternatively, FVRPs arise whenever some elements of a given problem involve uncertainty, subjectivity, ambiguity and vagueness (Teodorovic & Pavkovic, 1996). For instance, in OVRPFDs, information about each customer’s demand is often not adequately precise. For example, based on experience, a customer’s demand can be approximated as ‘‘around 50 units,” ‘‘between 20 and 60 units,” and so on. Generally, we can use fuzzy variables to deal with these uncertain parameters. In fact, Teodorovic and Pavkovic (1996)used fuzzy variables with regard to VRPs, and Cheng and Gen (1995) used a genetic algorithm to solve VRPs with fuzzy due time. Finally, moreover, Lai, Liu, and Peng (2003) modeled VRPs with fuzzy travel times using fuzzy programming with a possibility measure; they then adopted a genetic algorithm to solve the model. Zheng and Liu (2006)researched VRPs with fuzzy travel times and presented a chance-constrained program (CCP) model with a credibility measure; they then integrated fuzzy simulation and a genetic algorithm (GA) to design a hybrid intelligent algorithm to solve the model. Note that FVRPs differ from their deterministic counterparts in several fundamental respects. The concept of a solution is different, as several fundamental properties of deterministic VRPs no longer hold in FVRPs. Thus, solution methodologies are considerably more intricate. To the best of our knowledge, no research has considered non-deterministic information in an OVRP framework. In this paper, we consider the situation in which the demands of customers are fuzzy variables in OVRP. We model the OVRPFD by using fuzzy programming with a credibility measure and adopt a hybrid differential algorithm to solve the model. Differential evolution (DE) is a novel evolutionary technique that was originally developed for continuous optimization (Storn, 1996; Storn & Price, 1995). DE is a population-based, global evolutionary algorithm, which uses both a simple operator to create new candidate solutions and a one-to-one competition scheme to select new candidates greedily. Due to its simple structure, easy implementation, quick convergence, and robustness, DE has been one of the best evolutionary algorithms for solving continuous problems in a variety of fields. Nevertheless, due to the fact that DE is continuous, the research on DE for combinatorial optimization is very limited. Obviously, it is difficult to apply DE to problems other than the continuous optimization on which its inventors originally focused. Recently, some researchers have used DE to address machine layout problems (Nearchou, 2006; Storn, 1996) and to solve manufacturing problems with mixed-integer, discrete variables (Lampinen & Zelinka, 1999). But to the best of our knowledge, there is no work on VRPs that uses differential evolution. In this paper, we will first adopt the differential evolution algorithm to solve routing problems; the proposed differential evolution algorithm also can solve other deterministic aspects of the OVRPFD insofar as it considers more complicated constraints. This paper is organized as follows: in Section 2, we explain some basic concepts regarding fuzzy theory. In Section 3, we introduce open vehicle routing problems with fuzzy demands and present a chance-constrained program (CCP) model with which we measure fuzzy events with credibility. Then in Section 4, we integrate both stochastic simulation and a differential evolution algo-
rithm to design a hybrid intelligent algorithm to solve the CCP model. In Section 5, we discuss two experiments in order to demonstrate the effectiveness of the hybrid intelligent algorithm. In the final section, we summarize the contributions of this paper. 2. Fuzzy credibility measure theory The concept of a fuzzy set was first discussed by Zadeh (1965) with regard to membership functions. Since then, the concept of fuzzy sets has been well developed and applied in a wide variety of real problems. In order to measure a fuzzy event, the term fuzzy variable was introduced by Kaufman (1975), while possibility measure theory of fuzzy variable was proposed by Zadeh (1978). Liu (2004) recently developed credibility theory. In this section, we introduce briefly some basic concepts in fuzzy measure theory. First, we will introduce the axioms of possibility measure theory (Liu, 2004), which form the basis of credibility measure theory. Let H be a nonempty set, and let P be the power set of H. Each element in P is called an event, and u is an empty set. In order to present an axiomatic definition of possibility, it is necessary to assign a number Pos{A} to each event A, which indicates the possibility that A will occur. In order to ensure that the number Pos{A} has certain mathematical properties that we would intuitively expect a possibility to have, we accept the following four axioms: Axiom 2.1. Pos{H} = 1 (Liu, 2004); Axiom 2.2. Pos{u} = 0; Axiom 2.3. For each Ai p(H), Pos{[iAi} = supiPos{Ai}; Axiom2.4. If Hi is a non-empty set, and the set function Posi{}, i = 1,2, . . . , n, satisfies the above three axioms, and H = H1 H2 sup Pos1{h1}^ Hn, then for each A e p(H), PosfAg ¼ ðh1 ;h2 ;...;hn Þ2A Pos2{h2}^. . .Posn{hn}. The above four axioms form the basis of credibility measure theory; thus, we can derive all concepts of credibility theory from them (Liu, 2004). Definition 2.5. Let (H, P(H), Pos) be a possibility space, and let A be a set in p(H). Then the necessity measure of A is defined by Nec{A} = 1 Pos{A}(Liu, 2004). Definition 2.6. Let (H, P(H), Pos) be a possibility space, and A be a set in p(H), then the credibility measure of A is defined by CrfAg ¼ 12 ðPosfAg þ NecfAgÞ (Liu, 2004). Obviously, a fuzzy event may not hold even though its possibility approaches 1 and such an event may hold even though its necessity is 0. However, a fuzzy event must hold if its credibility is 1, and it must fail if its credibility is 0 (Liu, 2004). The credibility measure is self-dual, and in the theory of fuzzy subsets, the law of credibility plays a role similar to that played by the law of probability in measurement theory for ordinary sets. Now let us consider a triangular fuzzy variable d = (d1, d2, d3) as the demand of a given customer such that d is described by its left boundary d1 and its right boundary d3. Thus, a dispatcher or analyst studying such a problem can subjectively estimate, based on his/ her experience, intuition and/or available data, that a customer’s demand will not be less than d1 or greater than d3. The value of d2 corresponds to a grade of membership of 1, which can also be determined by a subjective estimate. Based on these definitions of possibility, necessity and credibility, we can derive:
Posfd rg ¼
8 > < 1; > :
d3 r d3 d2
0;
if ; if if
r d2 d2 r d3 r d3
ð1Þ
2407
E. Cao, M. Lai / Expert Systems with Applications 37 (2010) 2405–2411
Necfd rg ¼
Crfd rg ¼
8 > < 1;
d2 r d2 d1
if ; if
r d1 d1 r d2
> : 0; if r d2 8 1; if r d1 > > > > 2d2 d1 r ; if d r d < 1 2 2ðd2 d1 Þ d3 r > ; > 2ðd d Þ > > : 3 2 0;
if
d2 r d3
if
r d3
ð2Þ
ð3Þ
3. The fuzzy chance-constrained program model of the OVRPFD We assume that: (a) each vehicle has a container with a physical limitation so that the total loading of each vehicle can not exceed its capacity C; (b) each vehicle has maximum distance constraints so that the total distance traveled by each vehicle can not exceed L; (c) a vehicle will be assigned for only one route on which there may be more than one customer; (d) a customer will be visited by one and only one vehicle; (e) each route begins at the distribution depot (0) and ends at the last customer; (f) the demand of each customer is a triangular fuzzy number d = (d1, d2, d3); (g) the distance between customer i and customer j is cij; (h) and there are k vehicle in the distribution depot. For convenience, the capacity of all vehicles is the same, and the maximum distance that each vehicle can travel is also equal; these constraints can be denoted respectively as C and L. We assume that all vehicles are filled with goods required by customers when they leave the distribution depot, and each vehicle must fulfill the prestated services and end at the last customer. We will denote a triangular fuzzy number representing demand at the ith customer by di = (d1i, d2i, d3i). After serving the first k customers, the demand a P vehicle is still able to fulfill is Q k ¼ C ki¼1 di . Qk is also a triangular fuzzy and thus requires the use of fuzzy arithmetic, where P P P Q k ¼ ðC ki¼1 d3i ; C ki¼1 d2i ; C ki¼1 d1i Þ ¼ ðq1;k ; q2;k ; q3;k Þ. Thus, we can evaluate the credibility that the demand of the next customers does not exceed a vehicle’s supply. That is,
Cr ¼ C r fdkþ1 Q k g ¼ C r fðd1;kþ1 q3;k ; d2;kþ1 q2;k ; d3;kþ1 q1;k ; Þ 0g 8 0; if d1;kþ1 q3;k > > > q3;k d1;kþ1 > < ; if d1;kþ1 q3;k ; d2;kþ1 q2;k 2ðq3;k d1;kþ1 þd2;kþ1 q2;k Þ ¼ d3;kþ1 q1;k 2ðd2;kþ1 q2;k Þ > if d2;kþ1 q2;k ; d3;kþ1 q1;k > þd q Þ 2ðq d > > : 2;k 2;kþ1 3;kþ1 1;k 1 if d3;kþ1 q1;k
ð4Þ
As we now know, if a vehicle’s remaining supply is relatively large and the demand at the next customer is relatively low, then the vehicle’s ‘‘chance” of being able to finish the next customer’s service should relatively favorable. That is to say, the greater the difference between the vehicle’s available supply and the demand at the next customer, the greater is our preference to send the vehicle to serve the next customer. We will describe the preference index by Cr, which denotes the magnitude of our preference to send the vehicle to the next customer after it has served the current customer according to formulation (4). Note that Cr e [0, 1]. When Cr = 0, we declare that the vehicle should terminate service at the current customer and return the depot. When Cr = 1, we are completely sure that the vehicle should serve the next customer. Let the dispatcher preference index equal Cr*, where Cr* e [0, 1]. Note that Cr* expresses the dispatcher’s attitude toward risk. When the dispatcher is not risk-averse, he/she will choose lower values of parameter Cr*, which indicates the dispatcher’s preference to use vehicles as much as possible, even though there is an increase in the number of cases in which the vehicle arrives at the next customer and is not able to carry out planned service due to lack of goods. Alternatively, if the dispatcher chooses a greater Cr* in order to maximize the occurrences in which the vehicle supply is greater
than the demand of future customers, then the dispatcher is considered risk-averse. So, according to the dispatcher preference index value and the credibility that the next customer’s demand does not exceed the vehicle supply, a decision must be made regarding whether to send the vehicle to the next customer or to return it to the depot to renew its supply of goods. In this paper, this decision is made as follows: if Cr P Cr* holds, then the vehicle is sent to the next customer; otherwise, the vehicle is returned to the depot to renew its supply. It then continues to service the next customer. We do not terminate the above process until the demands of all customers are fulfilled. Vehicle routes are designed in advance by applying this algorithm. However, the actual value of demand at a customer is only known when the vehicle reaches that customer. Due to the uncertainty of demand, a vehicle might not be able to service a customer once it arrives due to insufficient supply with regard to the planned route. In this paper, the vehicle is assumed to return to the depot to pick up goods and then returns to the next customer whose demands remain unfulfilled. The vehicle then continues to service the rest of the planned route; accordingly, the vehicle travels additional distance due to this ‘‘failure” en route. Thus, when evaluating the planned route, we must consider the additional distance that the vehicle may travel due to ‘‘failure” along the route. Parameter Cr*, which is subjectively determined in advance, has an extremely great impact on both the final distance of the planned routes and on the additional distance covered by vehicles due to ‘‘failures.” As already mentioned, lower values of parameter Cr* express the dispatcher’s desire to use vehicles’ remaining supplies the best that he/she can. These values result in shorter planned route length. But lower values of parameter Cr* increase the number of situations in which vehicles arrive at a customer and are unable to service it due to small available supply, thereby increasing the total distance they cover due to these ‘‘failures.” We use stochastic simulation to evaluate the additional distance due to failures in following section. Alternatively, higher values of parameter Cr* are characterized by lower utilization of remaining vehicle supply along the planned routes and less additional distance as a result of fewer failures. Thus, the problem that logically arises is deciding which value of parameter Cr* will result in the shortest route, including the planned route as well as additional distance due to failures. We assume the decision variable xijk = 1 if arc(i, j) belongs to the route operated by vehicle k, and otherwise xijk = 0; yik = 1 if customer i is serviced by vehicle k, otherwise yik = 0. The corresponding chance-constrained problem (CCP), that is, the mathematical formulation of OVRPFD based on credibility theory, is as follows:
min
k X n X n X k¼1 i¼0
cij xijk
min c0 s:t: Cr
n X
ð5Þ
j¼0
ð6Þ
! di yik C
Cr ;
k ¼ 1; . . . ; k
ð7Þ
i¼1 k X
yik ¼ 1;
i ¼ 1; 2; . . . ; n;
ð8Þ
k¼1 n X i¼0 n X
xijk
n X
xjik ¼ 0;
j ¼ 0; 1; . . . ; n; k ¼ 1; . . . ; k;
ð9Þ
i¼0
xj0k 1;
k ¼ 1; 2; . . . ; k
ð10Þ
j¼1 n X i¼0
xijk ¼ yjk ;
j ¼ 1; . . . ; n; k ¼ 1; . . . ; k
ð11Þ
2408
E. Cao, M. Lai / Expert Systems with Applications 37 (2010) 2405–2411
n X n X i¼0
cij xijk L;
k ¼ 1; . . . ; k
ð12Þ
j¼0
xijk 2 f0; 1g; yik 2 f0; 1g;
i; j ¼ 0; 1; . . . ; n; k ¼ 1; . . . ; k
ð13Þ
Function (5) minimizes total planned travel distance. Function (6) minimizes total additional travel distance due to route failures; note that c0 can be obtain using the stochastic simulation algorithm in Section 4.1. The sum of planned and additional distance is obtained using the improved differential evolution algorithm in Section 4.2. Chance constraint (7) assures that all customers are visited given a vehicle’s capacity within a certain confidence level. Constraint (8) ensures that each customer is visited by exactly one vehicle. Constraint (9) guarantees that the same vehicle arrives and departs from each customer it serves. Restriction (10) requires that at most k vehicles are used. Restriction (11) requires that each vehicle starts at the distribution depot and ends at the last customer, and it also indicates the relation between two decision variables. Restriction (12) involves the maximum distance constraints; L is the upper limit on the total distances transported by a vehicle in any given section of the route. Finally, constraint (13) defines the nature of the decision variable. 4. A hybrid intelligent algorithm In this paper, we design a hybrid intelligent algorithm integrating stochastic simulation and a differential evolution algorithm to solve the above fuzzy chance-constrained problem. For a given value of dispatcher preference index C r , we adopt a triangular fuzzy number within vehicle capacity to represent demand at each customer; the real value of demand at a given customer is a real number within fuzzy boundaries that is made known when the vehicle reaches the customer. First, we apply stochastic simulation to simulate the real value of demands at each customer in order to derive the planned distance as well as the additional distance due to ‘‘route failures”. Then, we apply an improved differential evolution algorithm to obtain the total expected distance for all vehicles, which is the total sum of the length of all planned routes plus additional distance covered by vehicles. The objective is to obtain the value of parameter Cr* that will result in the least total sum. 4.1. Stochastic simulation algorithm Because each customer’s demand is a triangular fuzzy number, we cannot deal with it directly as with a deterministic number. The real value of a customer’s demand when the vehicle reaches that customer can be considered, however, to be a deterministic number using simulation. For each feasible planned route, we derive an approximate estimate about additional distances (c0 ) due to route failures using a stochastic simulation algorithm. We summarize the stochastic simulation as follows. Algorithm 1. Step1: For each customer, estimate the additional distances by simulating ‘‘actual” demands. The ‘‘actual” demands are generated through the following processes: (1) randomly generate a real number x in the interval between the left and right boundaries of the triangular fuzzy number, which represents customer demand, and compute its membership u; (2) generate a random number a, a e [0, 1]; (3) compare a with u; if a 6 u, then the ‘‘actual” customer demand is set equal to x, while in the opposite case, if a > u, customer does not equal x. In this case, random numbers x and a are generated again
and again until x and a are found that satisfy a 6 u; (4) repeat (1)–(3), and terminate the process when each customer has a simulated ‘‘actual” demand quantity. Step2: Move vehicles along the routes designed using credibility theory and deposit goods at each customer. Calculate the additional distance due to route failures in terms of ‘‘actual” demand. Step3: Repeat Step 1 and Step 2 M times. Step4: Compute the average value of additional distance based on M simulations, and regard it as the additional distance. 4.2. Differential evolution algorithm The differential evolution (DE) algorithm grew out of Price’s attempts to solve the Chebychev Polynomial Fitting Problem that had been posed to him by Storn (1996). DE combines simple arithmetic operators with classical operators such as crossover, mutation and selection to evolve from a randomly generated starting population to a final solution. The crucial idea behind DE involves the new scheme for generating trial parameter vectors. DE generates new parameter vectors by adding the weighted difference vector between two population members to a third member. If the resulting vector yields a lower objective function value than a predetermined population member, the newly-generated vector replaces the vector with which it was compared. In addition, the best parameter vector is evaluated for every generation in order to keep track of the progress that is made during the minimization process. DE uses a rather greedy and less stochastic approach in problem-solving as compared to evolution algorithms. These simple, yet powerful and straightforward features make it very attractive in numerical optimization. When applying DE, we face many difficulties. Key difficulties include dealing with fuzzy numbers regarding customer demand and extending the canonical DE to optimize integer variable problems. The customer permutation-based encoding scheme (Dantzing & Ramser, 1959; Cheng & Gen, 1995; Zheng & Liu, 2006) has been widely used for VRP. However, due to the continuous nature of DE, the standard encoding scheme of DE cannot be directly adopted for VRP. So, the important objective in applying DE to VRP is finding a suitable map between customer sequence and individuals, both of which are continuous vectors in DE. The improved differential evolution algorithm is briefly described below insofar as it is used in order to effectively solve OVRPFDs. More detailed strategies regarding DE can be founded in DE (2008) and Storn and Price (1995). 4.2.1. Coding and initiation population As with all evolutionary optimization algorithms, DE is a novel parallel direct search method that works with a population of solutions, rather than with a single solution, for the optimization problem. Like conventional GA for the VRP (Gen & Cheng, 1997; Zheng &Liu, 2006), a chromosome Chrom(i,:) is a sequence (or permutation) of n customers. In adopting this ordinal number encoding method, we note that the initial population is chosen randomly and that population size is NP. NP does not change during the minimization process, and each individual is an n-dimensional solution vector that is a permutation possible customer numbers, i.e.,
Chromði; :Þ ¼ randpermðnÞ;
i ¼ 1; 2; . . . ; NP;
where n is the number of customers. We check the capacity chance constraint (7) and distance constraint (12) at the first customer (or gene) of the permutation (or chromosome). If the constraints are not violated, the next customer is considered. If constraints are violated at a certain customer, we use another vehicle starting with this customer and repeat the above process until all customers
2409
E. Cao, M. Lai / Expert Systems with Applications 37 (2010) 2405–2411
are serviced. That is to say, the demand of each customer must line up with the capacity chance constraint in a route, and in any route, the total traveled distance can not exceed the maximum distance of the combined distances of each vehicle, For instance, for 10 customers, a randomly generated permutation is 1 3 6 8 9 5 4 10 2 7,1 3 6 8 9 5 4 10 2 7, which can be interpreted as r = 3 feasible routes: 0-1-36-0, 0-8-9-5-0, and 0-4-10-2-7-0. If k >¼ r, then this chromosome is legal; otherwise, it is illegal. In order to prevent an illegal chromosome from entering the next generation with great probability, a penalty function is designed. R is the total distance traveled for a given chromosome, and let m ¼ r k. If r > k, then m > 0, and R = R+Z m, where Z is a very large integer. If r < k, then m = 0. The fitness function can be expressed as f = 1/(R+Z m). We summarize the above process as follows: Algorithm 2 Step1: Randomly generate chromosomes of customer arrangements. Step2: Choose the first customer of a chromosome. According to customer demand and remaining vehicle supply, we can derive Cr from formula (4). For a dispatcher preference index value Cr*, if Cr P Cr*, then the customer is assigned to the current vehicle. Otherwise, another vehicle is used to service this customer. Step3: Delete the first customer from the chromosome. Step4: Repeat Step 2 and Step 3. If all customers have been assigned to routes, we have generated a feasible chromosome. Step5: Repeat Steps 1–4 for a given population size NP. 4.2.2. Mutation operation Mutation generates chromosome of offspring based on parent gene difference. It is an operation that adds a vector differential to a population vector of individuals according to the following equation:
v ði; :Þ ¼ Chromðc; :Þ þ F ½Chromða; :Þ Chromðb; :Þ;
ð14Þ
where a, b, c, i e [1, NP] are randomly selected and are mutually different from each other and also different from the running index i, i.e., a – b – c – i. These values refer to three randomly-chosen population vectors. The scaling factor F e [0, 2] is a real constant, and it controls the amplification of the differential variation Chrom(a,:) Chrom(b,:). Because we adopt an ordinal encoding scheme and because each chromosome represents a permutation of the customers, each gene stands for a customer. When the offspring gene lies outside a range allowable for formula (14), we consider an auxiliary operator based on integer order criterion (IOR), similar to the largest-order-value (LOV) rules that Qian, Wang, Huang, et al. (2009) proposed. The largest gene of offspring gives the largest customer ordinal number n, the second is evaluated as n 1, and the rest may be deduced by analogy. We can prove that this operator equates to an affine transform (Onwubolu & Davendra, 2006). For example, an offspring chromosome is [7.1, 1.3, 5.6, 2.5, 3.7, 0, 3.3, 5.4]. The number of customers is 8, and we obtain an offspring chromosome [1, 5, 2, 6, 3, 4, 7, 8] using IOR. We adopt the ordinal encoding method and an auxiliary operator based on IOR to improve the mutation operation, which is a crucial reason why differential evolution can solve vehicle routing problems. If other more common rounded number methods would have been adopted, we could not have solved vehicle routing problems. 4.2.3. Crossover operation In order to increase the potential diversity of the perturbed parameter vectors, a crossover operation is introduced after the
mutation operation. The crossover operation is employed to generate a temporary or trial vector by replacing certain parameters of the target vector with corresponding parameters of a randomly generated donor vector. The trial vector is formed using the following equation:
( trialði; jÞGþ1 ¼
v ði; jÞGþ1 ;
randj CR or j ¼ randnðiÞ
Chromði; jÞG ; otherwise
;
ð15Þ
where randj is a random value within the interval [0, 1]; and randn(i) is a randomly chosen index from the set of customers. We obtain the ordinal valued vector using the mutation operation, and therefore all vectors generated from formula (15) are integer-valued vectors. Note that the index randn(i) refers to a randomly chosen vector parameter, and it is used to ensure that at least one vector parameter of each individual trial vector trial(i,:)G+1 differs from its counterpart in the previous generation Chrom(i,:)G. In another words, a certain sequence of vector elements of trial(i,:)G+1 is identical to the elements of v(i,:)G+1 if a stochastic number smaller than CR is generated; otherwise, the other elements of trial(i,:)G+1 acquire the original values of Chrom(i,:)G. In formula (15), G is the number of current iterations, and CR e [0, 1] is the crossover probability factor. In order to improve the population’s diversity and its ability to break away from the local optimum, we present a new self-adapting crossover probability factor. The crossover probability (CR) is time-varying, as it changes from small to large at G iterations:
CR ¼ CRmin þ G
CRmax CRmin ; MAXGEN
ð16Þ
where CRmin is the proposed minimum crossover probability; CRmax is the maximum crossover probability; and MAXGEN is the number of maximum iteration. In the early stages of evolution, crossover probability is smaller, which can improve global searching capabilities; in the later stages of evolution, crossover probability is larger, which can improve local searching capabilities. 4.2.4. Estimation and selection operations According to formula (5), we derive the total planned length of all routes (c) and the additional distance covered due to route failures (c0 ) by using stochastic simulation in Section 4.1. We define the fitness value as the reciprocal of the total length of all routes, 1 i.e., f ¼ 1R ¼ ðcþc 0 . The DE selection scheme differs from other evolutionary algorithms. The fitness of the temporary individual trial(i,:)G+1 is compared with that of its target individual Chrom(i,:)G in each generation. The one with higher fitness value will propagate the population of the next generation. The population of the next generation Chrom(i,:)G+1 (i = 1,2, ... , NP) is obtained by using the following greedy selection criterion:
( Chromði; : Þ
Gþ1
¼
trialði; : ÞGþ1 ;
f ðChromði; : ÞG Þ < f ðtrialði; : ÞGþ1 Þ
Chromði; :Þ; otherwise ð17Þ
Usually, the main DE parameters are fewer than those in other evolutionary algorithms. Note that the performance of a DE algorithm depends on three parameters: the population size NP, the scaling factor F and the crossover probability CR. Practical advice on how to select control parameters NP, F and CR can be found in Storn (1996) and Storn and Price (1995).
5. Numerical experiments Now we will offer some examples to illustrate the models that we have just discussed and show how the hybrid intelligent algorithm works. Two types of experimental conditions are created
2410
E. Cao, M. Lai / Expert Systems with Applications 37 (2010) 2405–2411
Table 1 The relative parameters by model and algorithm. n
C
30
8
5000
k
G
M
L
NP
CRmin
CRmax
F
30
100
100
2000
60
0.3
0.9
0.5
Planned distances Additional distances Total distances
4500 4000 3500 3000
based on the size of the problem (in this case, the number of customers). We assume that there are 30 customers and one depot for the small-size problem and 100 customers and one depot for the large-size problem. In each experiment, the coordinates of all customers and the depot are randomly generated in [100 100], and the fuzzy demands of customers are randomly generated so that they are triangular fuzzy numbers within vehicle capacity C. We estimate the additional distances due to route failures using the stochastic simulation algorithm in Section 4.1, and we derive the planned distances and total distances using the improved differential evolution algorithm in Section 4.2. The relative parameters for the small-size problem are listed in Table 1. For the 100-customer problem, n = 100, k ¼ 100, NP = 100, and the other parameters are identical to those in Table 1. We also derive planned distances, additional distances and total distances, and we reveal the influence of the dispatcher preference index Cr* on these distances. The hybrid intelligent algorithm was encoded using MATLAB 7.0. The value of dispatcher preference index Cr* varied within the interval [0, 1] with a step of 0.1. The average computational results for 10 iterations are given in Tables 2 and 3 for the small-size and large-size problems, respectively. Figs. 1 and 2, respectively, show tendencies regarding planned distances, additional distances due to failures, and total distances that vehicles covered as the dispatcher preference index varied. In Tables 2 and 3, as the dispatcher preference index Cr* rose, the total distances of planned routes strictly increased, while additional distance due to failures strictly decreased. When Cr* 6 0.6,
2500 2000 1500 1000 500 0
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
*
Fig. 1. Changes in distance as Cr varies when n = 30.
18000 Planned distances Additional distances Total distances
16000 14000 12000 10000 8000 6000 4000 2000
Table 2 Average results with different Cr* when n = 30.
0
Cr*
Planned distances (103)
Additional distances (103)
Total distances (103)
0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0
1.4652 2.2873 2.4603 2.5384 2.7646 2.8647 3.0063 3.1991 3.2966 3.5203 3.6982
3.2409 1.4562 1.0969 0.7709 0.4981 0.2733 0.1217 0.0512 0.0008 0.0000 0
4.7061 3.7436 3.5572 3.3093 3.2627 3.1380 3.1280 3.2503 3.2974 3.5203 3.6982
Table 3 Average results with different Cr* when n = 100. Cr*
Planned distances (104)
Additional distances (104)
Total distances (104)
0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0
0.5322 0.7336 0.7899 0.8482 0.8851 0.9516 1.0094 1.0499 1.0985 1.1177 1.1388
1.0899 0.5373 0.3947 0.2792 0.1834 0.1204 0.0396 0.0115 0.0007 0 0
1.6221 1.2709 1.1845 1.1275 1.0685 1.0720 1.0491 1.0614 1.0992 1.1177 1.1388
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
Fig. 2. Changes in distances as Cr* varies when n = 100.
the increase in planned distances is smaller than the decrease in the additional distances, and so the total distance is strictly decreasing as Cr* increases from 0 to 0.6. However, when Cr* > 0.6, of the increase in planned distances is strictly larger than the decrease in additional distances, and so the total distance is almost strictly increasing as Cr* increases from 0.6 to 1. When the dispatcher preference value equals 0.6, the total distance covered by a given vehicle is minimized. As a consequence, lower values of Cr* model a subjective desire to use vehicle capacity the best that we can. These values correspond to routes with shorter planned distances. However, lower values of Cr* increase the number of cases in which a vehicle is unable to service customer demands, thereby increasing the total additional distance they cover due to ‘‘failures.” Higher values of Cr* are characterized by lower utilization of vehicle capacity along the planned routes and less additional distance due to fewer failures. Therefore, we conclude that the dispatcher preference index should be approximately 0.6. 6. Conclusion This paper contributed to understanding vehicle routing problem with fuzzy demands in the following respects: (a) a
E. Cao, M. Lai / Expert Systems with Applications 37 (2010) 2405–2411
chance-constrained model of OVRPFD was proposed based on credibility theory; (b) stochastic simulation and differential evolution algorithm were integrated to design a hybrid intelligent algorithm used to solve this problem by focusing on minimizing total distance; (c) the dispatcher preference index greatly influenced the length of planned routes and the additional distances covered by vehicles due to failures. Note that a ‘‘best” value of parameter Cr* was derived using the proposed hybrid algorithm; and (d) the effectiveness of the hybrid intelligent algorithm was demonstrated through numerical examples. In the further, we may consider other non-deterministic side constraints in OVRP such as time constraints, including traveling time, as well as customer locations. We may also adopt other meta-heuristics to solve these problems, such as the tabu search, the ant colony algorithm, and particle swarm optimization. Acknowledgment This research was supported by the Specialized Research Fund for Doctoral Program of Higher Education of China under Grant 20050532029. References Bertsimas, D. J. (1992). A vehicle routing problem with stochastic demand. Operational Research, 40, 574–585. Brandão, J. (2004). A tabu search heuristic algorithm for open vehicle routing problem. European Journal of Operational Research, 157, 552–564. Cheng, R., & Gen, M. (1995). Vehicle routing problem with fuzzy due-time using genetic algorithm. Japanese Journal of Fuzzy Theory and Systems, 7(5), 1050–1061. Dantzing, G., & Ramser, J. (1959). The truck dispatching problem. Management Science, 10(6), 80–91. DE Home page. (2008). http://www.icsi.berkeley.edu/~storn/code.html. Dror, M., Laporte, G., & Trudeau, P. (1989). Vehicle routing with stochastic demands: properties and solution frameworks. Transportation Science, 23, 166–176. Fu, Z., Eglese, R., & Li, L. Y. O. (2005). A new tabu search heuristic for the open vehicle routing problem. Journal of the Operational Research Society, 56, 267–274. Gen, M., & Cheng, R. (1997). Genetic algorithm and engineering design. New York: John Wiley and Sons, Inc..
2411
Gendreau, M., Laporate, G., & Seguin, R. (1996). Stochastic vehicle routing. European Journal of Operational Research, 88, 3–12. Kaufman, A. (1975). Introduction to the theory of fuzzy subsets. New York: Academic Press. Lai, K. K., Liu, B., & Peng, J. (2003). Vehicle routing problem with fuzzy travel times and its genetic algorithm. Technical Report. Lampinen, J., & Zelinka, I. (1999). Mechanical engineering design optimization by differential evolution. In D. Corne, M. Dorigo, & F. Glover (Eds.), New ideas in optimization (pp. 127–146). London, UK: McGraw-Hill. Li, F., Golden, B., & Wasil, E. (2007). The open vehicle routing problem: algorithms, large-scale test problems, and computational results. Computers and Operations Research, 38(10), 2918–2930. Liu, B. (2004). Uncertain theory: An introduce to its axiomatic foundations. Berlin: Springer. Liu, B., & Lai, K. K. (2002). Stochastic programming models for vehicle routing problems. Asian Information Science Life, 1(1), 13–28. Nearchou, A. C. (2006). Meta-heuristics from nature for the loop layout design problem. International Journal of Production Economics, 101(2), 312–328. Onwubolu, G., & Davendra, D. (2006). Scheduling flow shops using differential evolution algorithm. European Journal of Operational Research, 171, 674–692. Pisinger, D., & Ropke, S. (2007). A general heuristic for vehicle routing problems. Computers and Operations Research, 38(4), 2403–2435. Qian, B., Wang, Ling, Huang, De-xian, et al. (2009). An effective hybrid DE-based algorithm for multi-objective flow shop scheduling with limited buffers. Computers and Operations Research, 36(1), 209–233. Reinelt, G. (1991). The traveling salesman: Computational solutions for TSP applications. Berlin: Springer-Verlag. Schrage, L. (1981). Formulation and structure of more complex/realistic routing and scheduling problems. Networks, 11, 229–232. Storn, R. (1996). Differential evolution design of an IIR-filter. In Proceedings IEEE conference evolutionary computation (pp. 268–273). Nagoya, Japan: IEEE. Storn, R., & Price, K. (1995). Differential evolution – a simple and efficient adaptive scheme for global optimization over continuous spaces. Technical Report TR95-012, ICSI, March, ftp.icsi.berkeley.edu. Tarantilis, C. D., Ioannou, G., Kiranoudis, C. T., & Prastacos, G. P. (2005). Solving the open vehicle routing problem via single parameter meta-heuristic algorithm. Journal of the Operational Research Society, 56, 588–596. Tarantilis, C. D., & Kiranoudis, C. T. (2002). Distribution of fresh meat. Journal of Food Engineering, 51, 85–91. Teodorovic, D., & Pavkovic, G. (1996). The fuzzy set theory approach to the vehicle routing problem when demand at nodes is uncertain. Fuzzy Sets and Systems, 82, 307–317. Zadeh, L. A. (1978). Fuzzy sets as a basis for a theory of possibility. Fuzzy Sets and Systems, 1, 3–28. Zadeh, L. A. (1965). Fuzzy sets. Information and Control, 8, 338–353. Zheng, Y., & Liu, B. (2006). Fuzzy vehicle routing model with credibility measure and its hybrid intelligent algorithm. Applied Mathematics and Computation, 176, 673–683.
Contents lists available at ScienceDirect
Expert Systems with Applications journal homepage: www.elsevier.com/locate/eswa
The open vehicle routing problem with fuzzy demands Cao Erbao *, Lai Mingyong College of Economics and Trade, Hunan University, Changsha 410079, China
a r t i c l e
i n f o
Keywords: Logistics distribution Open vehicle routing problem Fuzzy credibility Stochastic simulation Differential evolution algorithm Optimization
a b s t r a c t According to the open vehicle routing problem (OVRP), a vehicle is not required to return to the distribution depot after servicing the last customer on its route. In this paper, the open vehicle routing problem with fuzzy demands (OVRPFD) is considered. A fuzzy chance-constrained program model is designed based on fuzzy credibility theory. Stochastic simulation and an improved differential evolution algorithm are integrated so as to use a hybrid intelligent algorithm to solve the OVRPFD model. The influence of the decision-maker’s preference on the final outcome of the problem is analyzed using stochastic simulation, and the range of possible preferences is calculated. Crown Copyright Ó 2009 Published by Elsevier Ltd. All rights reserved.
1. Introduction The open vehicle routing problem (OVRP) analyzes efficient routes with minimum total cost for a fleet of vehicles that serve some commodity to a given number of customers. Each customer is visited exactly once by one vehicle, while vehicle activity is bounded by capacity constraints, duration constraints and time constraints. Either each route is a sequence of customers that starts at the depot and finishes at one of the customers to whom goods are delivered, or each route is a sequence of customers that begins at a certain customer and ends at the distribution depot, where goods are gathered. The OVRP differs from the well-known vehicle routing problem (VRP) (Dantzing & Ramser, 1959) in that the vehicles do not necessarily return to their original locations after delivering goods to customers; if they do, they must visit the same customers in the reverse order. The OVRP is encountered in practice in many contexts, such as home delivery of packages and newspapers. Contractors who are not employees of the delivery company use their own vehicles and do not return to the depot; the ‘‘missing” node on a route can be driver’s home or a parking lot where the vehicle stays overnight. The major difference in theory between the OVRP and the VRP is that the routes in the OVRP consist of Hamiltonian paths originating at the depot and terminating at one of the customers, while the routes in the VRP are Hamiltonian cycles. In other words, the best Hamiltonian path is NP-hard, since the Hamiltonian path problem is equivalent to the traveling salesperson problem, which is known to be NP-hard (Reinelt, 1991). The best Hamiltonian path problem with a fixed source node must be solved for each vehicle in the OVRP, and OVRP solutions involve finding the best Hamilto* Corresponding author. Tel.: +86 7318646317; fax: +86 7318684825. E-mail address: [email protected] (E. Cao).
nian path for each set of customers assigned to a vehicle. Consequently, the OVRP is also a NP-hard problem. The OVRP has received sparse attention in the literature compared to the VRP. While the earliest description of OVRP offered by Schrage (1981) appeared in the literature over 20 years ago, OVRPs have just recently attracted the attention of practitioners and researchers. From the early 1980s to the late 1990s, the OVRP received very little attention in the operations research literature. However, since 2000, several researchers have used various heuristics and meta-heuristics to solve the OVRP with some success. Brandão (2004) and Fu, Eglese, and Li (2005) implemented a tabu search (TS) heuristic to solve the OVRP with constraints on vehicle capacity and maximum route length. Tarantilis and Kiranoudis (2002) considered a real-life fresh meat distribution problem, formulating it as a multi-depot OVRP; the problem was solved using a metaheuristic method called the list-based threshold-accepting (LBTA) algorithm. Tarantilis, Ioannou, Kiranoudis, and Prastacos (2005) solved the OVRP by adopting the LBTA algorithm they proposed for the solution of the multi-depot OVRP (Tarantilis & Kiranoudis, 2002). Pisinger and Ropke (2007) presented an adaptive largeneighborhood search heuristic to solve an OVRP. Li, Golden, and Wasil (2007) proposed a record-to-record travel heuristic and a deterministic variant of simulated annealing to solve the OVRP. However, traditional studies of the OVRP as well as the VRP have assumed that all information is deterministic, including customer information, vehicles information, roads information, and dispatcher information. As such, the proposed algorithms were only used to solve deterministic circumstances. That is, all the above-mentioned papers assumed that the demands of all customers visited on its route by any vehicle were deterministic. Actually, in some new systems, it is hard to describe the parameters of the VRP as deterministic, because there exists much
0957-4174/$ - see front matter Crown Copyright Ó 2009 Published by Elsevier Ltd. All rights reserved. doi:10.1016/j.eswa.2009.07.021
2406
E. Cao, M. Lai / Expert Systems with Applications 37 (2010) 2405–2411
uncertainty in data such as customer demands, traveling time and the set of customers to be visited. We call these problems nondeterministic VRPs. Given this aspect of VRPs, a consideration of stochastic vehicle routing problems (SVRP) and fuzzy vehicle routing problems (FVRP) may be useful. SVRPs arise whenever some elements of a given problem are random. Common examples include stochastic demands and stochastic travel times. Sometimes, the set of customers to be visited is not known with certainty. In this case, each customer has a certain probability of being visited. Researchers have developed many models and algorithms for SVRPs (Bertsimas, 1992; Dror, Laporte, & Trudeau, 1989; Gendreau, Laporate, & Seguin, 1996; Liu & Lai, 2002). Alternatively, FVRPs arise whenever some elements of a given problem involve uncertainty, subjectivity, ambiguity and vagueness (Teodorovic & Pavkovic, 1996). For instance, in OVRPFDs, information about each customer’s demand is often not adequately precise. For example, based on experience, a customer’s demand can be approximated as ‘‘around 50 units,” ‘‘between 20 and 60 units,” and so on. Generally, we can use fuzzy variables to deal with these uncertain parameters. In fact, Teodorovic and Pavkovic (1996)used fuzzy variables with regard to VRPs, and Cheng and Gen (1995) used a genetic algorithm to solve VRPs with fuzzy due time. Finally, moreover, Lai, Liu, and Peng (2003) modeled VRPs with fuzzy travel times using fuzzy programming with a possibility measure; they then adopted a genetic algorithm to solve the model. Zheng and Liu (2006)researched VRPs with fuzzy travel times and presented a chance-constrained program (CCP) model with a credibility measure; they then integrated fuzzy simulation and a genetic algorithm (GA) to design a hybrid intelligent algorithm to solve the model. Note that FVRPs differ from their deterministic counterparts in several fundamental respects. The concept of a solution is different, as several fundamental properties of deterministic VRPs no longer hold in FVRPs. Thus, solution methodologies are considerably more intricate. To the best of our knowledge, no research has considered non-deterministic information in an OVRP framework. In this paper, we consider the situation in which the demands of customers are fuzzy variables in OVRP. We model the OVRPFD by using fuzzy programming with a credibility measure and adopt a hybrid differential algorithm to solve the model. Differential evolution (DE) is a novel evolutionary technique that was originally developed for continuous optimization (Storn, 1996; Storn & Price, 1995). DE is a population-based, global evolutionary algorithm, which uses both a simple operator to create new candidate solutions and a one-to-one competition scheme to select new candidates greedily. Due to its simple structure, easy implementation, quick convergence, and robustness, DE has been one of the best evolutionary algorithms for solving continuous problems in a variety of fields. Nevertheless, due to the fact that DE is continuous, the research on DE for combinatorial optimization is very limited. Obviously, it is difficult to apply DE to problems other than the continuous optimization on which its inventors originally focused. Recently, some researchers have used DE to address machine layout problems (Nearchou, 2006; Storn, 1996) and to solve manufacturing problems with mixed-integer, discrete variables (Lampinen & Zelinka, 1999). But to the best of our knowledge, there is no work on VRPs that uses differential evolution. In this paper, we will first adopt the differential evolution algorithm to solve routing problems; the proposed differential evolution algorithm also can solve other deterministic aspects of the OVRPFD insofar as it considers more complicated constraints. This paper is organized as follows: in Section 2, we explain some basic concepts regarding fuzzy theory. In Section 3, we introduce open vehicle routing problems with fuzzy demands and present a chance-constrained program (CCP) model with which we measure fuzzy events with credibility. Then in Section 4, we integrate both stochastic simulation and a differential evolution algo-
rithm to design a hybrid intelligent algorithm to solve the CCP model. In Section 5, we discuss two experiments in order to demonstrate the effectiveness of the hybrid intelligent algorithm. In the final section, we summarize the contributions of this paper. 2. Fuzzy credibility measure theory The concept of a fuzzy set was first discussed by Zadeh (1965) with regard to membership functions. Since then, the concept of fuzzy sets has been well developed and applied in a wide variety of real problems. In order to measure a fuzzy event, the term fuzzy variable was introduced by Kaufman (1975), while possibility measure theory of fuzzy variable was proposed by Zadeh (1978). Liu (2004) recently developed credibility theory. In this section, we introduce briefly some basic concepts in fuzzy measure theory. First, we will introduce the axioms of possibility measure theory (Liu, 2004), which form the basis of credibility measure theory. Let H be a nonempty set, and let P be the power set of H. Each element in P is called an event, and u is an empty set. In order to present an axiomatic definition of possibility, it is necessary to assign a number Pos{A} to each event A, which indicates the possibility that A will occur. In order to ensure that the number Pos{A} has certain mathematical properties that we would intuitively expect a possibility to have, we accept the following four axioms: Axiom 2.1. Pos{H} = 1 (Liu, 2004); Axiom 2.2. Pos{u} = 0; Axiom 2.3. For each Ai p(H), Pos{[iAi} = supiPos{Ai}; Axiom2.4. If Hi is a non-empty set, and the set function Posi{}, i = 1,2, . . . , n, satisfies the above three axioms, and H = H1 H2 sup Pos1{h1}^ Hn, then for each A e p(H), PosfAg ¼ ðh1 ;h2 ;...;hn Þ2A Pos2{h2}^. . .Posn{hn}. The above four axioms form the basis of credibility measure theory; thus, we can derive all concepts of credibility theory from them (Liu, 2004). Definition 2.5. Let (H, P(H), Pos) be a possibility space, and let A be a set in p(H). Then the necessity measure of A is defined by Nec{A} = 1 Pos{A}(Liu, 2004). Definition 2.6. Let (H, P(H), Pos) be a possibility space, and A be a set in p(H), then the credibility measure of A is defined by CrfAg ¼ 12 ðPosfAg þ NecfAgÞ (Liu, 2004). Obviously, a fuzzy event may not hold even though its possibility approaches 1 and such an event may hold even though its necessity is 0. However, a fuzzy event must hold if its credibility is 1, and it must fail if its credibility is 0 (Liu, 2004). The credibility measure is self-dual, and in the theory of fuzzy subsets, the law of credibility plays a role similar to that played by the law of probability in measurement theory for ordinary sets. Now let us consider a triangular fuzzy variable d = (d1, d2, d3) as the demand of a given customer such that d is described by its left boundary d1 and its right boundary d3. Thus, a dispatcher or analyst studying such a problem can subjectively estimate, based on his/ her experience, intuition and/or available data, that a customer’s demand will not be less than d1 or greater than d3. The value of d2 corresponds to a grade of membership of 1, which can also be determined by a subjective estimate. Based on these definitions of possibility, necessity and credibility, we can derive:
Posfd rg ¼
8 > < 1; > :
d3 r d3 d2
0;
if ; if if
r d2 d2 r d3 r d3
ð1Þ
2407
E. Cao, M. Lai / Expert Systems with Applications 37 (2010) 2405–2411
Necfd rg ¼
Crfd rg ¼
8 > < 1;
d2 r d2 d1
if ; if
r d1 d1 r d2
> : 0; if r d2 8 1; if r d1 > > > > 2d2 d1 r ; if d r d < 1 2 2ðd2 d1 Þ d3 r > ; > 2ðd d Þ > > : 3 2 0;
if
d2 r d3
if
r d3
ð2Þ
ð3Þ
3. The fuzzy chance-constrained program model of the OVRPFD We assume that: (a) each vehicle has a container with a physical limitation so that the total loading of each vehicle can not exceed its capacity C; (b) each vehicle has maximum distance constraints so that the total distance traveled by each vehicle can not exceed L; (c) a vehicle will be assigned for only one route on which there may be more than one customer; (d) a customer will be visited by one and only one vehicle; (e) each route begins at the distribution depot (0) and ends at the last customer; (f) the demand of each customer is a triangular fuzzy number d = (d1, d2, d3); (g) the distance between customer i and customer j is cij; (h) and there are k vehicle in the distribution depot. For convenience, the capacity of all vehicles is the same, and the maximum distance that each vehicle can travel is also equal; these constraints can be denoted respectively as C and L. We assume that all vehicles are filled with goods required by customers when they leave the distribution depot, and each vehicle must fulfill the prestated services and end at the last customer. We will denote a triangular fuzzy number representing demand at the ith customer by di = (d1i, d2i, d3i). After serving the first k customers, the demand a P vehicle is still able to fulfill is Q k ¼ C ki¼1 di . Qk is also a triangular fuzzy and thus requires the use of fuzzy arithmetic, where P P P Q k ¼ ðC ki¼1 d3i ; C ki¼1 d2i ; C ki¼1 d1i Þ ¼ ðq1;k ; q2;k ; q3;k Þ. Thus, we can evaluate the credibility that the demand of the next customers does not exceed a vehicle’s supply. That is,
Cr ¼ C r fdkþ1 Q k g ¼ C r fðd1;kþ1 q3;k ; d2;kþ1 q2;k ; d3;kþ1 q1;k ; Þ 0g 8 0; if d1;kþ1 q3;k > > > q3;k d1;kþ1 > < ; if d1;kþ1 q3;k ; d2;kþ1 q2;k 2ðq3;k d1;kþ1 þd2;kþ1 q2;k Þ ¼ d3;kþ1 q1;k 2ðd2;kþ1 q2;k Þ > if d2;kþ1 q2;k ; d3;kþ1 q1;k > þd q Þ 2ðq d > > : 2;k 2;kþ1 3;kþ1 1;k 1 if d3;kþ1 q1;k
ð4Þ
As we now know, if a vehicle’s remaining supply is relatively large and the demand at the next customer is relatively low, then the vehicle’s ‘‘chance” of being able to finish the next customer’s service should relatively favorable. That is to say, the greater the difference between the vehicle’s available supply and the demand at the next customer, the greater is our preference to send the vehicle to serve the next customer. We will describe the preference index by Cr, which denotes the magnitude of our preference to send the vehicle to the next customer after it has served the current customer according to formulation (4). Note that Cr e [0, 1]. When Cr = 0, we declare that the vehicle should terminate service at the current customer and return the depot. When Cr = 1, we are completely sure that the vehicle should serve the next customer. Let the dispatcher preference index equal Cr*, where Cr* e [0, 1]. Note that Cr* expresses the dispatcher’s attitude toward risk. When the dispatcher is not risk-averse, he/she will choose lower values of parameter Cr*, which indicates the dispatcher’s preference to use vehicles as much as possible, even though there is an increase in the number of cases in which the vehicle arrives at the next customer and is not able to carry out planned service due to lack of goods. Alternatively, if the dispatcher chooses a greater Cr* in order to maximize the occurrences in which the vehicle supply is greater
than the demand of future customers, then the dispatcher is considered risk-averse. So, according to the dispatcher preference index value and the credibility that the next customer’s demand does not exceed the vehicle supply, a decision must be made regarding whether to send the vehicle to the next customer or to return it to the depot to renew its supply of goods. In this paper, this decision is made as follows: if Cr P Cr* holds, then the vehicle is sent to the next customer; otherwise, the vehicle is returned to the depot to renew its supply. It then continues to service the next customer. We do not terminate the above process until the demands of all customers are fulfilled. Vehicle routes are designed in advance by applying this algorithm. However, the actual value of demand at a customer is only known when the vehicle reaches that customer. Due to the uncertainty of demand, a vehicle might not be able to service a customer once it arrives due to insufficient supply with regard to the planned route. In this paper, the vehicle is assumed to return to the depot to pick up goods and then returns to the next customer whose demands remain unfulfilled. The vehicle then continues to service the rest of the planned route; accordingly, the vehicle travels additional distance due to this ‘‘failure” en route. Thus, when evaluating the planned route, we must consider the additional distance that the vehicle may travel due to ‘‘failure” along the route. Parameter Cr*, which is subjectively determined in advance, has an extremely great impact on both the final distance of the planned routes and on the additional distance covered by vehicles due to ‘‘failures.” As already mentioned, lower values of parameter Cr* express the dispatcher’s desire to use vehicles’ remaining supplies the best that he/she can. These values result in shorter planned route length. But lower values of parameter Cr* increase the number of situations in which vehicles arrive at a customer and are unable to service it due to small available supply, thereby increasing the total distance they cover due to these ‘‘failures.” We use stochastic simulation to evaluate the additional distance due to failures in following section. Alternatively, higher values of parameter Cr* are characterized by lower utilization of remaining vehicle supply along the planned routes and less additional distance as a result of fewer failures. Thus, the problem that logically arises is deciding which value of parameter Cr* will result in the shortest route, including the planned route as well as additional distance due to failures. We assume the decision variable xijk = 1 if arc(i, j) belongs to the route operated by vehicle k, and otherwise xijk = 0; yik = 1 if customer i is serviced by vehicle k, otherwise yik = 0. The corresponding chance-constrained problem (CCP), that is, the mathematical formulation of OVRPFD based on credibility theory, is as follows:
min
k X n X n X k¼1 i¼0
cij xijk
min c0 s:t: Cr
n X
ð5Þ
j¼0
ð6Þ
! di yik C
Cr ;
k ¼ 1; . . . ; k
ð7Þ
i¼1 k X
yik ¼ 1;
i ¼ 1; 2; . . . ; n;
ð8Þ
k¼1 n X i¼0 n X
xijk
n X
xjik ¼ 0;
j ¼ 0; 1; . . . ; n; k ¼ 1; . . . ; k;
ð9Þ
i¼0
xj0k 1;
k ¼ 1; 2; . . . ; k
ð10Þ
j¼1 n X i¼0
xijk ¼ yjk ;
j ¼ 1; . . . ; n; k ¼ 1; . . . ; k
ð11Þ
2408
E. Cao, M. Lai / Expert Systems with Applications 37 (2010) 2405–2411
n X n X i¼0
cij xijk L;
k ¼ 1; . . . ; k
ð12Þ
j¼0
xijk 2 f0; 1g; yik 2 f0; 1g;
i; j ¼ 0; 1; . . . ; n; k ¼ 1; . . . ; k
ð13Þ
Function (5) minimizes total planned travel distance. Function (6) minimizes total additional travel distance due to route failures; note that c0 can be obtain using the stochastic simulation algorithm in Section 4.1. The sum of planned and additional distance is obtained using the improved differential evolution algorithm in Section 4.2. Chance constraint (7) assures that all customers are visited given a vehicle’s capacity within a certain confidence level. Constraint (8) ensures that each customer is visited by exactly one vehicle. Constraint (9) guarantees that the same vehicle arrives and departs from each customer it serves. Restriction (10) requires that at most k vehicles are used. Restriction (11) requires that each vehicle starts at the distribution depot and ends at the last customer, and it also indicates the relation between two decision variables. Restriction (12) involves the maximum distance constraints; L is the upper limit on the total distances transported by a vehicle in any given section of the route. Finally, constraint (13) defines the nature of the decision variable. 4. A hybrid intelligent algorithm In this paper, we design a hybrid intelligent algorithm integrating stochastic simulation and a differential evolution algorithm to solve the above fuzzy chance-constrained problem. For a given value of dispatcher preference index C r , we adopt a triangular fuzzy number within vehicle capacity to represent demand at each customer; the real value of demand at a given customer is a real number within fuzzy boundaries that is made known when the vehicle reaches the customer. First, we apply stochastic simulation to simulate the real value of demands at each customer in order to derive the planned distance as well as the additional distance due to ‘‘route failures”. Then, we apply an improved differential evolution algorithm to obtain the total expected distance for all vehicles, which is the total sum of the length of all planned routes plus additional distance covered by vehicles. The objective is to obtain the value of parameter Cr* that will result in the least total sum. 4.1. Stochastic simulation algorithm Because each customer’s demand is a triangular fuzzy number, we cannot deal with it directly as with a deterministic number. The real value of a customer’s demand when the vehicle reaches that customer can be considered, however, to be a deterministic number using simulation. For each feasible planned route, we derive an approximate estimate about additional distances (c0 ) due to route failures using a stochastic simulation algorithm. We summarize the stochastic simulation as follows. Algorithm 1. Step1: For each customer, estimate the additional distances by simulating ‘‘actual” demands. The ‘‘actual” demands are generated through the following processes: (1) randomly generate a real number x in the interval between the left and right boundaries of the triangular fuzzy number, which represents customer demand, and compute its membership u; (2) generate a random number a, a e [0, 1]; (3) compare a with u; if a 6 u, then the ‘‘actual” customer demand is set equal to x, while in the opposite case, if a > u, customer does not equal x. In this case, random numbers x and a are generated again
and again until x and a are found that satisfy a 6 u; (4) repeat (1)–(3), and terminate the process when each customer has a simulated ‘‘actual” demand quantity. Step2: Move vehicles along the routes designed using credibility theory and deposit goods at each customer. Calculate the additional distance due to route failures in terms of ‘‘actual” demand. Step3: Repeat Step 1 and Step 2 M times. Step4: Compute the average value of additional distance based on M simulations, and regard it as the additional distance. 4.2. Differential evolution algorithm The differential evolution (DE) algorithm grew out of Price’s attempts to solve the Chebychev Polynomial Fitting Problem that had been posed to him by Storn (1996). DE combines simple arithmetic operators with classical operators such as crossover, mutation and selection to evolve from a randomly generated starting population to a final solution. The crucial idea behind DE involves the new scheme for generating trial parameter vectors. DE generates new parameter vectors by adding the weighted difference vector between two population members to a third member. If the resulting vector yields a lower objective function value than a predetermined population member, the newly-generated vector replaces the vector with which it was compared. In addition, the best parameter vector is evaluated for every generation in order to keep track of the progress that is made during the minimization process. DE uses a rather greedy and less stochastic approach in problem-solving as compared to evolution algorithms. These simple, yet powerful and straightforward features make it very attractive in numerical optimization. When applying DE, we face many difficulties. Key difficulties include dealing with fuzzy numbers regarding customer demand and extending the canonical DE to optimize integer variable problems. The customer permutation-based encoding scheme (Dantzing & Ramser, 1959; Cheng & Gen, 1995; Zheng & Liu, 2006) has been widely used for VRP. However, due to the continuous nature of DE, the standard encoding scheme of DE cannot be directly adopted for VRP. So, the important objective in applying DE to VRP is finding a suitable map between customer sequence and individuals, both of which are continuous vectors in DE. The improved differential evolution algorithm is briefly described below insofar as it is used in order to effectively solve OVRPFDs. More detailed strategies regarding DE can be founded in DE (2008) and Storn and Price (1995). 4.2.1. Coding and initiation population As with all evolutionary optimization algorithms, DE is a novel parallel direct search method that works with a population of solutions, rather than with a single solution, for the optimization problem. Like conventional GA for the VRP (Gen & Cheng, 1997; Zheng &Liu, 2006), a chromosome Chrom(i,:) is a sequence (or permutation) of n customers. In adopting this ordinal number encoding method, we note that the initial population is chosen randomly and that population size is NP. NP does not change during the minimization process, and each individual is an n-dimensional solution vector that is a permutation possible customer numbers, i.e.,
Chromði; :Þ ¼ randpermðnÞ;
i ¼ 1; 2; . . . ; NP;
where n is the number of customers. We check the capacity chance constraint (7) and distance constraint (12) at the first customer (or gene) of the permutation (or chromosome). If the constraints are not violated, the next customer is considered. If constraints are violated at a certain customer, we use another vehicle starting with this customer and repeat the above process until all customers
2409
E. Cao, M. Lai / Expert Systems with Applications 37 (2010) 2405–2411
are serviced. That is to say, the demand of each customer must line up with the capacity chance constraint in a route, and in any route, the total traveled distance can not exceed the maximum distance of the combined distances of each vehicle, For instance, for 10 customers, a randomly generated permutation is 1 3 6 8 9 5 4 10 2 7,1 3 6 8 9 5 4 10 2 7, which can be interpreted as r = 3 feasible routes: 0-1-36-0, 0-8-9-5-0, and 0-4-10-2-7-0. If k >¼ r, then this chromosome is legal; otherwise, it is illegal. In order to prevent an illegal chromosome from entering the next generation with great probability, a penalty function is designed. R is the total distance traveled for a given chromosome, and let m ¼ r k. If r > k, then m > 0, and R = R+Z m, where Z is a very large integer. If r < k, then m = 0. The fitness function can be expressed as f = 1/(R+Z m). We summarize the above process as follows: Algorithm 2 Step1: Randomly generate chromosomes of customer arrangements. Step2: Choose the first customer of a chromosome. According to customer demand and remaining vehicle supply, we can derive Cr from formula (4). For a dispatcher preference index value Cr*, if Cr P Cr*, then the customer is assigned to the current vehicle. Otherwise, another vehicle is used to service this customer. Step3: Delete the first customer from the chromosome. Step4: Repeat Step 2 and Step 3. If all customers have been assigned to routes, we have generated a feasible chromosome. Step5: Repeat Steps 1–4 for a given population size NP. 4.2.2. Mutation operation Mutation generates chromosome of offspring based on parent gene difference. It is an operation that adds a vector differential to a population vector of individuals according to the following equation:
v ði; :Þ ¼ Chromðc; :Þ þ F ½Chromða; :Þ Chromðb; :Þ;
ð14Þ
where a, b, c, i e [1, NP] are randomly selected and are mutually different from each other and also different from the running index i, i.e., a – b – c – i. These values refer to three randomly-chosen population vectors. The scaling factor F e [0, 2] is a real constant, and it controls the amplification of the differential variation Chrom(a,:) Chrom(b,:). Because we adopt an ordinal encoding scheme and because each chromosome represents a permutation of the customers, each gene stands for a customer. When the offspring gene lies outside a range allowable for formula (14), we consider an auxiliary operator based on integer order criterion (IOR), similar to the largest-order-value (LOV) rules that Qian, Wang, Huang, et al. (2009) proposed. The largest gene of offspring gives the largest customer ordinal number n, the second is evaluated as n 1, and the rest may be deduced by analogy. We can prove that this operator equates to an affine transform (Onwubolu & Davendra, 2006). For example, an offspring chromosome is [7.1, 1.3, 5.6, 2.5, 3.7, 0, 3.3, 5.4]. The number of customers is 8, and we obtain an offspring chromosome [1, 5, 2, 6, 3, 4, 7, 8] using IOR. We adopt the ordinal encoding method and an auxiliary operator based on IOR to improve the mutation operation, which is a crucial reason why differential evolution can solve vehicle routing problems. If other more common rounded number methods would have been adopted, we could not have solved vehicle routing problems. 4.2.3. Crossover operation In order to increase the potential diversity of the perturbed parameter vectors, a crossover operation is introduced after the
mutation operation. The crossover operation is employed to generate a temporary or trial vector by replacing certain parameters of the target vector with corresponding parameters of a randomly generated donor vector. The trial vector is formed using the following equation:
( trialði; jÞGþ1 ¼
v ði; jÞGþ1 ;
randj CR or j ¼ randnðiÞ
Chromði; jÞG ; otherwise
;
ð15Þ
where randj is a random value within the interval [0, 1]; and randn(i) is a randomly chosen index from the set of customers. We obtain the ordinal valued vector using the mutation operation, and therefore all vectors generated from formula (15) are integer-valued vectors. Note that the index randn(i) refers to a randomly chosen vector parameter, and it is used to ensure that at least one vector parameter of each individual trial vector trial(i,:)G+1 differs from its counterpart in the previous generation Chrom(i,:)G. In another words, a certain sequence of vector elements of trial(i,:)G+1 is identical to the elements of v(i,:)G+1 if a stochastic number smaller than CR is generated; otherwise, the other elements of trial(i,:)G+1 acquire the original values of Chrom(i,:)G. In formula (15), G is the number of current iterations, and CR e [0, 1] is the crossover probability factor. In order to improve the population’s diversity and its ability to break away from the local optimum, we present a new self-adapting crossover probability factor. The crossover probability (CR) is time-varying, as it changes from small to large at G iterations:
CR ¼ CRmin þ G
CRmax CRmin ; MAXGEN
ð16Þ
where CRmin is the proposed minimum crossover probability; CRmax is the maximum crossover probability; and MAXGEN is the number of maximum iteration. In the early stages of evolution, crossover probability is smaller, which can improve global searching capabilities; in the later stages of evolution, crossover probability is larger, which can improve local searching capabilities. 4.2.4. Estimation and selection operations According to formula (5), we derive the total planned length of all routes (c) and the additional distance covered due to route failures (c0 ) by using stochastic simulation in Section 4.1. We define the fitness value as the reciprocal of the total length of all routes, 1 i.e., f ¼ 1R ¼ ðcþc 0 . The DE selection scheme differs from other evolutionary algorithms. The fitness of the temporary individual trial(i,:)G+1 is compared with that of its target individual Chrom(i,:)G in each generation. The one with higher fitness value will propagate the population of the next generation. The population of the next generation Chrom(i,:)G+1 (i = 1,2, ... , NP) is obtained by using the following greedy selection criterion:
( Chromði; : Þ
Gþ1
¼
trialði; : ÞGþ1 ;
f ðChromði; : ÞG Þ < f ðtrialði; : ÞGþ1 Þ
Chromði; :Þ; otherwise ð17Þ
Usually, the main DE parameters are fewer than those in other evolutionary algorithms. Note that the performance of a DE algorithm depends on three parameters: the population size NP, the scaling factor F and the crossover probability CR. Practical advice on how to select control parameters NP, F and CR can be found in Storn (1996) and Storn and Price (1995).
5. Numerical experiments Now we will offer some examples to illustrate the models that we have just discussed and show how the hybrid intelligent algorithm works. Two types of experimental conditions are created
2410
E. Cao, M. Lai / Expert Systems with Applications 37 (2010) 2405–2411
Table 1 The relative parameters by model and algorithm. n
C
30
8
5000
k
G
M
L
NP
CRmin
CRmax
F
30
100
100
2000
60
0.3
0.9
0.5
Planned distances Additional distances Total distances
4500 4000 3500 3000
based on the size of the problem (in this case, the number of customers). We assume that there are 30 customers and one depot for the small-size problem and 100 customers and one depot for the large-size problem. In each experiment, the coordinates of all customers and the depot are randomly generated in [100 100], and the fuzzy demands of customers are randomly generated so that they are triangular fuzzy numbers within vehicle capacity C. We estimate the additional distances due to route failures using the stochastic simulation algorithm in Section 4.1, and we derive the planned distances and total distances using the improved differential evolution algorithm in Section 4.2. The relative parameters for the small-size problem are listed in Table 1. For the 100-customer problem, n = 100, k ¼ 100, NP = 100, and the other parameters are identical to those in Table 1. We also derive planned distances, additional distances and total distances, and we reveal the influence of the dispatcher preference index Cr* on these distances. The hybrid intelligent algorithm was encoded using MATLAB 7.0. The value of dispatcher preference index Cr* varied within the interval [0, 1] with a step of 0.1. The average computational results for 10 iterations are given in Tables 2 and 3 for the small-size and large-size problems, respectively. Figs. 1 and 2, respectively, show tendencies regarding planned distances, additional distances due to failures, and total distances that vehicles covered as the dispatcher preference index varied. In Tables 2 and 3, as the dispatcher preference index Cr* rose, the total distances of planned routes strictly increased, while additional distance due to failures strictly decreased. When Cr* 6 0.6,
2500 2000 1500 1000 500 0
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
*
Fig. 1. Changes in distance as Cr varies when n = 30.
18000 Planned distances Additional distances Total distances
16000 14000 12000 10000 8000 6000 4000 2000
Table 2 Average results with different Cr* when n = 30.
0
Cr*
Planned distances (103)
Additional distances (103)
Total distances (103)
0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0
1.4652 2.2873 2.4603 2.5384 2.7646 2.8647 3.0063 3.1991 3.2966 3.5203 3.6982
3.2409 1.4562 1.0969 0.7709 0.4981 0.2733 0.1217 0.0512 0.0008 0.0000 0
4.7061 3.7436 3.5572 3.3093 3.2627 3.1380 3.1280 3.2503 3.2974 3.5203 3.6982
Table 3 Average results with different Cr* when n = 100. Cr*
Planned distances (104)
Additional distances (104)
Total distances (104)
0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0
0.5322 0.7336 0.7899 0.8482 0.8851 0.9516 1.0094 1.0499 1.0985 1.1177 1.1388
1.0899 0.5373 0.3947 0.2792 0.1834 0.1204 0.0396 0.0115 0.0007 0 0
1.6221 1.2709 1.1845 1.1275 1.0685 1.0720 1.0491 1.0614 1.0992 1.1177 1.1388
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
Fig. 2. Changes in distances as Cr* varies when n = 100.
the increase in planned distances is smaller than the decrease in the additional distances, and so the total distance is strictly decreasing as Cr* increases from 0 to 0.6. However, when Cr* > 0.6, of the increase in planned distances is strictly larger than the decrease in additional distances, and so the total distance is almost strictly increasing as Cr* increases from 0.6 to 1. When the dispatcher preference value equals 0.6, the total distance covered by a given vehicle is minimized. As a consequence, lower values of Cr* model a subjective desire to use vehicle capacity the best that we can. These values correspond to routes with shorter planned distances. However, lower values of Cr* increase the number of cases in which a vehicle is unable to service customer demands, thereby increasing the total additional distance they cover due to ‘‘failures.” Higher values of Cr* are characterized by lower utilization of vehicle capacity along the planned routes and less additional distance due to fewer failures. Therefore, we conclude that the dispatcher preference index should be approximately 0.6. 6. Conclusion This paper contributed to understanding vehicle routing problem with fuzzy demands in the following respects: (a) a
E. Cao, M. Lai / Expert Systems with Applications 37 (2010) 2405–2411
chance-constrained model of OVRPFD was proposed based on credibility theory; (b) stochastic simulation and differential evolution algorithm were integrated to design a hybrid intelligent algorithm used to solve this problem by focusing on minimizing total distance; (c) the dispatcher preference index greatly influenced the length of planned routes and the additional distances covered by vehicles due to failures. Note that a ‘‘best” value of parameter Cr* was derived using the proposed hybrid algorithm; and (d) the effectiveness of the hybrid intelligent algorithm was demonstrated through numerical examples. In the further, we may consider other non-deterministic side constraints in OVRP such as time constraints, including traveling time, as well as customer locations. We may also adopt other meta-heuristics to solve these problems, such as the tabu search, the ant colony algorithm, and particle swarm optimization. Acknowledgment This research was supported by the Specialized Research Fund for Doctoral Program of Higher Education of China under Grant 20050532029. References Bertsimas, D. J. (1992). A vehicle routing problem with stochastic demand. Operational Research, 40, 574–585. Brandão, J. (2004). A tabu search heuristic algorithm for open vehicle routing problem. European Journal of Operational Research, 157, 552–564. Cheng, R., & Gen, M. (1995). Vehicle routing problem with fuzzy due-time using genetic algorithm. Japanese Journal of Fuzzy Theory and Systems, 7(5), 1050–1061. Dantzing, G., & Ramser, J. (1959). The truck dispatching problem. Management Science, 10(6), 80–91. DE Home page. (2008). http://www.icsi.berkeley.edu/~storn/code.html. Dror, M., Laporte, G., & Trudeau, P. (1989). Vehicle routing with stochastic demands: properties and solution frameworks. Transportation Science, 23, 166–176. Fu, Z., Eglese, R., & Li, L. Y. O. (2005). A new tabu search heuristic for the open vehicle routing problem. Journal of the Operational Research Society, 56, 267–274. Gen, M., & Cheng, R. (1997). Genetic algorithm and engineering design. New York: John Wiley and Sons, Inc..
2411
Gendreau, M., Laporate, G., & Seguin, R. (1996). Stochastic vehicle routing. European Journal of Operational Research, 88, 3–12. Kaufman, A. (1975). Introduction to the theory of fuzzy subsets. New York: Academic Press. Lai, K. K., Liu, B., & Peng, J. (2003). Vehicle routing problem with fuzzy travel times and its genetic algorithm. Technical Report. Lampinen, J., & Zelinka, I. (1999). Mechanical engineering design optimization by differential evolution. In D. Corne, M. Dorigo, & F. Glover (Eds.), New ideas in optimization (pp. 127–146). London, UK: McGraw-Hill. Li, F., Golden, B., & Wasil, E. (2007). The open vehicle routing problem: algorithms, large-scale test problems, and computational results. Computers and Operations Research, 38(10), 2918–2930. Liu, B. (2004). Uncertain theory: An introduce to its axiomatic foundations. Berlin: Springer. Liu, B., & Lai, K. K. (2002). Stochastic programming models for vehicle routing problems. Asian Information Science Life, 1(1), 13–28. Nearchou, A. C. (2006). Meta-heuristics from nature for the loop layout design problem. International Journal of Production Economics, 101(2), 312–328. Onwubolu, G., & Davendra, D. (2006). Scheduling flow shops using differential evolution algorithm. European Journal of Operational Research, 171, 674–692. Pisinger, D., & Ropke, S. (2007). A general heuristic for vehicle routing problems. Computers and Operations Research, 38(4), 2403–2435. Qian, B., Wang, Ling, Huang, De-xian, et al. (2009). An effective hybrid DE-based algorithm for multi-objective flow shop scheduling with limited buffers. Computers and Operations Research, 36(1), 209–233. Reinelt, G. (1991). The traveling salesman: Computational solutions for TSP applications. Berlin: Springer-Verlag. Schrage, L. (1981). Formulation and structure of more complex/realistic routing and scheduling problems. Networks, 11, 229–232. Storn, R. (1996). Differential evolution design of an IIR-filter. In Proceedings IEEE conference evolutionary computation (pp. 268–273). Nagoya, Japan: IEEE. Storn, R., & Price, K. (1995). Differential evolution – a simple and efficient adaptive scheme for global optimization over continuous spaces. Technical Report TR95-012, ICSI, March, ftp.icsi.berkeley.edu. Tarantilis, C. D., Ioannou, G., Kiranoudis, C. T., & Prastacos, G. P. (2005). Solving the open vehicle routing problem via single parameter meta-heuristic algorithm. Journal of the Operational Research Society, 56, 588–596. Tarantilis, C. D., & Kiranoudis, C. T. (2002). Distribution of fresh meat. Journal of Food Engineering, 51, 85–91. Teodorovic, D., & Pavkovic, G. (1996). The fuzzy set theory approach to the vehicle routing problem when demand at nodes is uncertain. Fuzzy Sets and Systems, 82, 307–317. Zadeh, L. A. (1978). Fuzzy sets as a basis for a theory of possibility. Fuzzy Sets and Systems, 1, 3–28. Zadeh, L. A. (1965). Fuzzy sets. Information and Control, 8, 338–353. Zheng, Y., & Liu, B. (2006). Fuzzy vehicle routing model with credibility measure and its hybrid intelligent algorithm. Applied Mathematics and Computation, 176, 673–683.