A genetic algorithm based approach to vehicle routing problem with simultaneous pick-up and deliveries

Computers & Industrial Engineering 62 (2012) 755–761

Contents lists available at SciVerse ScienceDirect

Computers & Industrial Engineering journal homepage: www.elsevier.com/locate/caie

A genetic algorithm based approach to vehicle routing problem with simultaneous pick-up and deliveries A. Serdar Tasan a,⇑, Mitsuo Gen b,c a

Department of Industrial Engineering, Dokuz Eylul University, Izmir, Turkey Fuzzy Logic Systems Institute (FLSI), Iizuka, Fukuoka, Japan c Hanyang University, Ansan, South Korea b

a r t i c l e

i n f o

Article history: Available online 29 November 2011 Keywords: Genetic algorithms Reverse logistics Vehicle routing problem with simultaneous pick-up and deliveries

a b s t r a c t The vehicle routing problem with simultaneous pick-up and deliveries, which considers simultaneous distribution and collection of goods to/from customers, is an extension of the capacitated vehicle routing problem. There are various real cases, where fleet of vehicles originated in a depot serves customers with pick-up and deliveries from/to their locations. Increasing importance of reverse logistics activities make it necessary to determine efficient and effective vehicle routes for simultaneous pick-up and delivery activities. The vehicle routing problem with simultaneous pick-up and deliveries is also NP-hard as a capacitated vehicle routing problem and this study proposes a genetic algorithm based approach to this problem. Computational example is presented with parameter settings in order to illustrate the proposed approach. Moreover, performance of the proposed approach is evaluated by solving several test problems. Ó 2011 Elsevier Ltd. All rights reserved.

1. Introduction The vehicle routing problem (VRP) was first introduced by Dantzig and Ramser (1959). The problem is concerned with delivering goods to a set of customers with known demands through vehicle routes that begin and finish at the depot with minimum cost. Capacitated vehicle routing problem (CVRP) is the most elementary version of the VRP and seeks a number of routes for m number of vehicles with Q units homogenous capacity to minimize total transportation cost of routes while satisfying the delivery demands of n number of customer nodes. Each route has to start and finish at the depot and each customer has to be visited exactly once by one vehicle. Details on VRP and its variants with formulation, and solution methods can be found in Toth and Vigo (2001). There are various real cases, where customers require simultaneous pick-up of goods from their location in addition to delivery of goods to their location and a fleet of vehicles originated in a depot serves customers with pick-up and deliveries from/to their locations. Increasing importance of reverse logistics activities make it necessary to determine efficient and effective vehicle routes for simultaneous pick-up and delivery activities. The vehicle routing problem with simultaneous pick-up and deliveries (VRP-SPD),

⇑ Corresponding author. Address: Dokuz Eylul Universitesi Muhendislik Fakultesi, Endustri Muh. Bolumu Tinaztepe Kampus Buca, Izmir 35160, Turkey. Tel.: +90 232 3017619; fax: +90 232 3017608. E-mail addresses: [email protected] (A.S. Tasan), gen@flsi.or.jp (M. Gen). 0360-8352/$ - see front matter Ó 2011 Elsevier Ltd. All rights reserved. doi:10.1016/j.cie.2011.11.025

which considers simultaneous distribution and collection of goods to/from customers, is an extension of the CVRP. In a broader perspective, the VRP can be considered as a combination of the bin packing problem (BPP) and the travelling salesman problem (TSP). Both BPP and TSP are NP-hard problems (Falkenauer, 1996; Wolsey, 1998); thus VRP and its variant VRPSPD are also NP-hard. During the last few decades, the importance of reverse logistics has shown a rapid increase due to the environmental and economical issues as well as legal obligations. One of the most critical issues that affect the performance of the reverse logistics is the routing of vehicles. However, efficient and effective vehicle routing requires high computational efforts. Increasing importance of reverse logistics activities and the computational complexity of VRP-SPD motivates this study. The purpose of this study is to introduce an alternative approach to solve VRP-SPD. More specifically, this study proposes a genetic algorithm (GA) based approach, which uses permutation based representation and ensures feasibility, for solving VRP-SPD more efficiently and effectively. In order to evaluate the performance of the developed GA based approach, several test problems from the literature are formulated as mixed integer linear programming (MILP) models. Results of the MILP solver are compared with the solution provided by the developed approach. The rest of the paper is organized as follows; in Section 2, mathematical formulation and details of VRP-SPD are presented. In Section 3, the details of the proposed approach is introduced. An

756

A.S. Tasan, M. Gen / Computers & Industrial Engineering 62 (2012) 755–761

delivery. VRP-B can be thought as a two separate CVRP problems; CVRP for the delivery customers and CVRP for the pickup (backhaul) customers. There are several studies related to VRP-B in the literature (e.g. Liu & Chung, 2009; Ropke & Pisinger, 2006). The simultaneous delivery and pick-up activities of VRP-SPD make this problem more difficult than VRP-B to solve. Additionally, feasibility issue is more complex in VRP-SPD. Recently, Gajpal and Abad (2009) applied ant colony optimization to VRP-SPD. Ai and Kachitvichyanukul (2009) proposed a particle swarm optimization for VRP-SPD with random key based representation. Zachariadis, Tarantilis, and Kiranoudis (2010) developed an adaptive memory methodology for this routing problem. In the following sections, VRP-SPD is formulated as a MILP model based on model given by Dethloff (2001). Fig. 1. Illustration of VRP-SPD. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

2.1. Notations Sets:

illustrative example is given in Section 4. The performance of the proposed approach is evaluated with several test problems and the results are discussed in Section 5. Finally, Section 6 concludes the paper with future directions.

J: Set of customer nodes. J0: Set of all nodes including depot. V: Set of vehicles. Parameters:

2. Vehicle routing problem with simultaneous pick-up and deliveries In VRP-SPD customers require simultaneous pick-up of goods from their location in addition to delivery of goods to their location as illustrated in Fig. 1. In Fig. 1, red dots represent customer nodes while blue rectangle represents the depot; red arrows represent deliveries; green arrows represent pick-ups in customer nodes; solid lines are the arcs between nodes. A number of routes are sought for v number of vehicles with C units homogenous capacity, that minimize the total transportation cost of routes while satisfying the pick-up and delivery demands of the m number of customer nodes simultaneously. Each vehicle has to depart and return to the depot and each customer node has to be visited exactly once by one vehicle. VRP-SPD was introduced into the literature by Min (1989). His study considered distribution problem of a public library with 22 branch libraries and 2 vehicles based on the central library. Vehicles serve branch libraries everyday with deliveries and pick-ups at each branch library. Vehicle routes were determined by a solution approach that based on clustering customers according to their demands and vehicle capacities first, and then solving TSP for each cluster. While determining routes for each cluster, an iterative procedure was used in order to satisfy the feasibility of the route. Dethloff (2001) emphasized the importance of VRP-SPD for reverse logistics activities and he investigated the relation between VRP-SPD and other types of routing problems. He proposed a mathematical formulation for this problem and developed an insertion based heuristic, which uses several insertion criteria. Nagy and Salhi (2005) also proposed insertion based heuristics in order to solve VRP-SPD. The basic steps of these heuristics are constructing partial routes for a set of customers, and then inserting the remaining customers to the existing route. They introduced the concepts of weak and strong feasibility and their proposed approach allowed infeasibilities to occur while searching towards strong feasibility. The proposed algorithm was extended to deal with problems with multiple depots. The vehicle routing problem with backhauls (VRP-B) is another extension of the VRP and considers both delivery and pick-ups like VRP-SPD. In the VRP-B, pick-up demands are met after finishing all

cap: Vehicle capacity. cij: Distance between nodes i and j i, j 2 J0, i – j. dj: Delivery amount demanded by customer node j 2 J. pj: Pick-up amount of customer node j 2 J. n: Number of nodes. M: Sufficiently large number (e.g. M can be calculated as the maximum value of either the total customer delivery and pick-up demands or the total distances between each nodes as given in (1)).

( ) X X X M ¼ max ðdj þ pj Þ; cij j2J

ð1Þ

i2J 0 j2J 0 ;i – j

Initialize Population

Decoding Routes

Selection

NO Genetic Operators NO

Survival

Terminate

YES END Fig. 2. Procedure of the proposed methodology.

757

A.S. Tasan, M. Gen / Computers & Industrial Engineering 62 (2012) 755–761

R1

R2

3. Proposed genetic algorithms based solution approach

R3

1 5 6 14 3 2 15 10 8 9 4 7 11 13 12 Fig. 3. Genetic representation of chromosome for 15 nodes.

Decision variables: 0

lv : Load of vth vehicle when leaving the depot. lj: Load of vehicle after served customer node j, j 2 J. sj: Variable used to avoid subtours, can be interpreted as position of node j 2 J in the route. xijv: Binary decision variable that indicates whether vth vehicle travels from i to j. 2.2. Mixed integer linear programming model for VRP-SPD The objective function of VRP-SPD model is given in (2) as the minimization of total travelled distance by vehicles those serve customer nodes.

XXX

Minimize

i2J 0 j2J 0

cij xijv

ð2Þ

v 2V

Subject to

XX

X

xijv ¼ 1;

xikv ¼

i2J 0

0

j2J

ð3Þ

v 2V

i2J 0

lv ¼

X

v 2V

xkjv ;

k 2 J;

dj xijv ;

v2V

ð4Þ

j2J 0

XX i2J 0

ð5Þ

j2J

0

lj P lv  dj þ pj  Mð1  x0jv Þ;

lj P li  dj þ pj  M 1 

X

j 2 J;

v2V

ð6Þ

! xijv ;

i 2 J;

j 2 J;

j–i

VRP-SPD belongs to the class of NP-hard problems, for that reason the exact solution methods become highly time-consuming as the problem instances increase in size. Therefore, due to the combinatorial nature of the VRP and the GAs’ efficiency in solving combinatorial problems, a GA based approach is developed to solve the vehicle routing problem. GAs can easily be adapted to various types of problems therefore many different GA approaches exist depending on the problems studied. There are several ways to maintain the population and several GA operators. However, all GA approaches must have a good genetic representation of the problem, an initial population generator, appropriate fitness function, and genetic operators such as crossover and mutation in order to work effectively. The generalized procedure of the GA approach in this study is shown in Fig. 2. The procedure for proposed GA based approach to solve VRPSPD (adapted from Gen, Cheng, and Lin (2008)) is: procedure: GA for VRP-SPD input: location data for nodes, delivery and pickup demands, vehicle capacity and GA parameters output: the best solution: (routes with shortest total travelled distance for serving customer nodes with simultaneous delivery and pick-ups) begin t 0; initialize P(t) by encoding routine; evaluate P(t) by decoding routine; while (not terminating condition) do create C(t) from P(t) by crossover routine; create C(t) from P(t) by mutation routine; evaluate C(t) by decoding routine; select P(t + 1) from P(t) and C(t) by selection routine; t t + 1; end output the best solution end

ð7Þ

v 2V

3.1. Genetic representation and encoding 0

lv 6 cap;

v2V

ð8Þ

lj 6 cap;

j2J

ð9Þ

sj P si þ 1  n 1 

X

! xijv ;

i 2 J;

j 2 J;

j–i

Instead of binary representation, permutation representation is used for genetic representation of the VRP-SPD as shown in Fig. 3. Direct representation is used to encode routing solution into chromosome. The length of the chromosome determined by the 100

ð10Þ

Customers

v 2V

sj P 0;

Depot

j2J

xijv 2 f0; 1g;

80

ð11Þ i 2 J0 ;

j 2 J0 ;

v 2V

60

ð12Þ

Constraint (3) ensures that each customer node is served exactly once. Constraint (4) guarantees that for each customer node, same vehicle arrives at and leaves this node. Initial vehicle loads are calculated as in (5), each vehicle’s initial load is the accumulated demand of all customer nodes assigned to this vehicle. (6) Balances the load of vehicles after vehicles visit the first customer node on their route. For other customer nodes the loads of vehicles are calculated as in (7) similarly through their routes. (8) and (9) Impose capacity constraint. Constraint (10) ensures sub tour elimination, and constraint (11) maintains non-negativity of sj and xijv is a binary decision variable as in (12).

Y 40

20

0

0

20

40

60

80

X Fig. 4. Locations of depot and customer nodes.

100

120

758

A.S. Tasan, M. Gen / Computers & Industrial Engineering 62 (2012) 755–761

D e m a n d

Customer Nodes Fig. 5. Delivery and pick-up demands of customer nodes.

3.3. Fitness value and decoding

Table 1 Computational results for the parameter settings. Popsize

Number of generations

I

II

III

IV

V

50 100 150

300 300 300

1296 1088 1023

1412 1276 1192

1480 1392 1308

2579 1917 2279

69.92 109.79 158.99

50 100 150

500 500 500

1240 1164 1079

1308 1230 1178

1374 1275 1294

3018 2341 2238

112.67 175.99 270.85

50 100 150

1000 1000 1000

1036 1101 1128

1167 1201 1168

1316 1259 1209

2502 2321 2345

240.15 357.42 545.35

number of customer nodes which are served by the vehicles. Routes are determined by the capacities of vehicles. 3.2. Initial population construction The desired number of individuals is generated to obtain a population with desired size. In the proposed methodology, initial population generation process is based on random permutation.

Due to the characteristic of VRP-SPD, the fitness value is calculated as shown in (2) that is the total distance travelled. The solutions are evaluated in terms of their fitness values which is identical to the fitness of individuals. The individuals with better fitness values survive while the ones with worse fitness values die. This means that the more costly solutions are removed from the population while others are remained (survived). Nagy and Salhi (2005) introduced the concepts of weak and strong feasibility in their study. In weakly feasible routes, neither the total pickup nor the total delivery load exceeds the vehicle capacity. In this study, weak feasibility is guaranteed with our decoding procedure and routes are determined according to the vehicle capacity, which ensure that total pickup and delivery loads through the route never exceeds the vehicle capacity. The load of vehicle, never exceeds its capacity after visiting each customer node in strongly feasible routes. To avoid strongly infeasible routes, a relatively high penalty cost is added to fitness value in the proposed approach. Penalty cost occurs each time when the load of vehicle exceeds capacity. In this study, while decoding procedure

I Best objective function value 3500

3000 2500

50-300

2000

100-300

1500 150-300

1000 IV Average objective function value of the population

500 0

II Average of the best values

50-500 100-500 150-500 50-1000 100-1000 150-1000

III Worst of the best values

Fig. 6. Results of the computations.

A.S. Tasan, M. Gen / Computers & Industrial Engineering 62 (2012) 755–761

120

100

N

Y 60 7

40

1 27

20

22

0

239 16

10

15 13

14 5

28

0

17 2012 29

34

24 4

R4

R5

25 21 8 33 31 6 18 3 26

80

changes will gradually add some new characteristics to the population, which could not be supplied by the crossover. In the study, partial-mapped crossover (PMX) and swap mutation (Gen et al., 2008) are used for genetic operations of permutation based chromosomes. After crossover and mutation, each offspring is evaluated in terms of fitness value mentioned before. The procedure of PMX is as follows (Gen et al., 2008):

R2

R1 R3

20

30

11 19

40

2

32

60

80

100

759

120

X Fig. 7. Routes for solution with total distance travelled 1023.

Procedure PMX input: Two parents output: Offspring step 1: Select two positions along the string uniformly random. The substring s defined by two positions are called mapping sections step 2: Exchange two substrings between parents to produce children step 3: Determine mapping relations between two mapping sections step 4: Legalize offspring with mapping relationship

guarantees weakly feasible routes, fitness value encourages strongly feasible routes. 4. Illustrative example and parameter settings 3.4. Selection

The main genetic operator is crossover, which simulates a reproduction between two parents. It works on a pair of solutions and recombines them in a certain way generating one or more offsprings. The offsprings share some of the characteristics of the parents and by this way the characteristic are passed onto the future generations. Crossover operator is not able to produce new characteristics. The other genetic operator is mutation, which is applied to a single solution with a certain probability. Mutation operator makes small random changes in the solution. These random

In this illustrative example, there are 34 customer nodes. These 34 customers are served with simultaneous pick-up and deliveries from/to depot by vehicles originated at depot. Vehicle capacity is 100 units. Location of depot and customer nodes and delivery and pick-up demands of customer nodes are shown in Figs. 4 and 5 respectively. Initially, to determine the appropriate population size and number of generations for GA, a pilot study is conducted. The proposed approach is applied with combinations of population size 50, 100, 150 and number of generations 300, 500 and 1000. The other parameters used in GA are crossover rate 0.80 and mutation rate 0.03. Results of several iterations are summarized in Table 1. For each combination, best objective function value (column I), average of the best values obtained through several runs (column II), worst of the best values obtained through several runs (column III), average objective function value of the population (column IV) as well as computation time as average CPU times (s) (column V) are given. Also Fig. 6 illustrates these results as a radar graph for column I to column IV. From these results, it can be concluded that population size 50 combined with number of generations 1000 results in good

Fig. 8. Best objective function value through the evolution process.

Fig. 9. Average objective function value through the evolution process.

The roulette wheel selection (RWS) is used in the proposed methodology, for probabilistically selection of individuals based on their performance. The roulette wheel selection method scales the fitness values of the members within the population so that the sum of the rescaled fitness values equals to 1. In the roulette wheel method, the probability of choosing an individual is directly proportional to its fitness value. 3.5. Genetic operators

760

A.S. Tasan, M. Gen / Computers & Industrial Engineering 62 (2012) 755–761

Table 2 The results of the computational experiments. Problem

Instance

Number of nodes

GA-Best

LB1-CVRP

LB2-CPLEX

LB = MAX(LB1, LB2)

UB(CPLEX)

CPLEX time (s)

GA time (s)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24

A-n32-k5.vrp A-n33-k5.vrp A-n33-k6.vrp A-n34-k5.vrp A-n36-k5.vrp A-n37-k5.vrp A-n37-k6.vrp A-n38-k5.vrp A-n39-k5.vrp A-n39-k6.vrp A-n45-k6.vrp A-n45-k7.vrp B-n31-k5.vrp B-n34-k5.vrp B-n35-k5.vrp B-n38-k6.vrp B-n39-k5.vrp B-n41-k6.vrp B-n43-k6.vrp B-n44-k7.vrp P-n16-k8.vrp P-n19-k2.vrp P-n20-k2.vrp P-n23-k8.vrp

32 33 33 34 36 37 37 38 39 39 45 45 31 34 35 38 39 41 43 44 16 19 20 23

957 781 798 923 1019 959 1115 974 1095 1257 1452 1397 736 865 1263 925 744 1105 1047 1126 479 240 222 569

784 661 742 778 799 669 949 730 822 831 944 1146 672 788 955 805 549 829 742 909 450 212 216 529

414.49 433.78 435.38 429.21 457.69 492.24 474.24 431.23 523.75 504.35 543.28 493.75 83.27 217.37 178.2 194.24 165.47 211.09 249.11 101.08 249.86 154.93 166.96 220.3

784 661 742 778 799 669 949 730 822 831 944 1146 672 788 955 805 549 829 742 909 450 212 216 529

1218.86 904.916 864.15 989.16 1246.58 992.55 1357.16 982.25 1192.65 1304.38 1627.91 1739.49 785.43 1004.02 1457.28 983.84 809.18 1111.01 1067 1227.69 478.63 257.18 266.14 611.98

1000 1000 1000 1000 1140 1000 1533 1000 1000 1047 1940 2660 1000 1000 1000 1000 1000 1000 2090 1000 1000 1000 1000 1000

226 233 266 264 275 287 287 280 286 302 326 333 241 264 274 284 296 311 315 332 181 159 159 235

solutions considering the average of best objective function values. Additionally combination of population size 150 with number of generation 300 is preferable when considering CPU time. The best solution found by the proposed approach is illustrated in Fig. 7. Total travelled distance by vehicles is 1023. Five vehicles served 34 customer nodes through: 1: 2: 3: 4: 5:

{depot, {depot, {depot, {depot, {depot,

5. Computational experiments and discussion In order to evaluate the performance of the proposed approach, twenty four medium-sized test problems are studied. These test problems are obtained by adding pick-up demand data to a selected set of CVRP problems from Augerat et al. (1995). The pick-up demands are generated as uniformly distributed values between [0, 26] for compatibility.

21, 8, 25, 31, 6, 18, 26, depot} 34, 4, 28, 5, 30, 11, depot} 12, 10, 16, 9, 23, 1, depot} 33, 3, 27, 24, 14, 22, 15, depot} 29, 17, 20, 7, 13, 32, 2, 19, depot}

2000 UB

GA

LB

1800 1600 1400

Objective Value

Route Route Route Route Route

For this example, best objective function and average objection function values through the evolution process are given in Figs. 8 and 9 respectively.

1200 1000 800 600 400 200 0 1

2

3

4

5

6

7

8

9

10 11 12 13 14 15 16 17 18 19 20 21 22 23 24

Problem # Fig. 10. Comparison of the results.

A.S. Tasan, M. Gen / Computers & Industrial Engineering 62 (2012) 755–761

During the computational experiments, parameter values which were previously determined in Section 4, are used as follows: Population size = 150. Maximum number of generations = 300. Crossover rate = 0.80. Mutation rate = 0.03. VRP-SPD problems are formulated as MILP model and CPLEX 10.0 solver (ILOG, 2006) is used for obtaining lower and upper bounds of the objective function. The best integer feasible solution given by CPLEX 10.0 after 1000 s are used as an upper bound (UB); in case of no feasible solution till 1000 s, the first feasible solution is considered as UB. CPLEX 10.0 produces lower bounds (LB-CPLEX) for the problems, however best known CVRP solutions given by Augerat et al. (1995) can be considered as alternative lower bounds (LB-CVRP). The results of the computations and CPU times with the details of test problems are shown in Table 2. Best known CVRP solutions is higher than lower bounds found by CPLEX for all test problems, so they are used as lower bounds LB = Max{LBCVRP, LB-CPLEX}. The LB, UB and results of GA based approach are shown in Fig. 10. As seen from Table 2 and Fig. 10, proposed GA provides solutions, which are between lower bound and upper bound for all of the 24 problems. In other words, solutions of the proposed GA based approach are better than the solutions provided by CPLEX 10.0 solver, which uses branch and bound algorithm for solving MILP models. In addition, the proposed approach gives better solutions within significantly shorter time frame when compared with CPLEX 10.0. Due to the NP-hard nature of VRP-SPD, MILP approach cannot solve the problem within an acceptable time frame. On the contrary, GAs perform well on NP-hard problems. Therefore, the proposed approach provided better solutions than MILP. Considering these results and CPU times, it can be stated that, GA based proposed approach is efficient and performs well. 6. Conclusion The VRP is a well known combinatorial optimization problem. Customers require simultaneous pick-up of goods from their location in addition to delivery of goods to their location in some cases. The VRP-SPD is an extension of the CVRP and considers simultaneous distribution and collection of goods to/from customers. A fleet of vehicles originated in a depot serves customers with pick-up and deliveries from/to their locations. VRP-SPD has been receiving growing attention due to the increasing importance of reverse logistics activities. Also VRPSPD has a combinatorial nature. Due to GA’s efficiency in solving complex optimization problems, a GA based approach for solving VRP-SPD was proposed in this study. Feasibility issue is more complex and hard to be ensured in VRP-SPD than in VRP-B due to the simultaneous delivery and pick-up activities as mentioned in Section 2. This study contributes to the VRP field by providing an efficient and effective GA based approach that produces highly feasible routes for VRP-SPD. Details of

761

the proposed approach were given in the previous sections, after introducing VRP-SPD and its mathematical formulation. Then, an illustrative example and parameter settings were presented followed by performance evaluation of the proposed approach with computational experiments. Twenty four CVRP test problems were transformed into VRP-SPD problems by adding pick-up demand data randomly to the selected problems. According to the results of computational experiments, it can be concluded that the proposed GA based approach both performs well and is efficient. The results of computational experiments prove that, the solutions of the proposed approach have both weak feasibility and strong feasibility. While decoding procedure guaranteed weak feasibility of solutions, fitness value enabled the proposed approach to produce strongly feasible routes. Direct genetic representation was used for encoding in this study. Further studies may use indirect representation methods such as random key based representation to improve solution quality and efficiency of the procedure. Additionally, the proposed approach may be applied to a real world routing problems with simultaneous pick-up and deliveries. Acknowledgment This research is partly supported by KOFST Project 2010-2011 at Hanyang University. References Ai, T. J., & Kachitvichyanukul, V. (2009). A particle swarm optimization for the vehicle routing problem with simultaneous pickup and delivery. Computers & Operations Research, 36, 1693–1702. Augerat, P., Belenguer, J. M., Benavent, E., Corberán, A., Naddef, D., & Rinaldi, G. (1995). Computational results with a branch and cut code for the capacitated vehicle routing problem. Research report 949-M. Grenoble, France: Universite Joseph Fourier. Dantzig, G. B., & Ramser, R. H. (1959). The truck dispatching problem. Management Science, 6, 80–91. Dethloff, J. (2001). Vehicle routing and reverse logistics: The vehicle routing problem with simultaneous delivery and pick-up. OR Spektrum, 23, 79–96. Falkenauer, E. (1996). A hybrid grouping genetic algorithm for bin packing. Journal of Heuristics, 2, 5–30. Gajpal, Y., & Abad, P. (2009). An ant colony system (ACS) for vehicle routing problem with simultaneous delivery and pickup. Computers & Operations Research, 36, 3215–3223. Gen, M., Cheng, R., & Lin, L. (2008). Network models and optimization: Multiobjective genetic algorithm approach. Springer-Verlag. ILOG (2006). ILOG CPLEX 10.0 user’s manual. ILOG SA and ILOG Inc. Liu, S. C., & Chung, C. H. (2009). A heuristic method for the vehicle routing problem with backhauls and inventory. Journal of Intelligent Manufacturing, 20, 29–42. Min, H. (1989). The multiple vehicle routing problem with simultaneous delivery and pick-up points. Transportation Research Part A: General, 23, 377–386. Nagy, G., & Salhi, S. (2005). Heuristic algorithms for single and multiple depot vehicle routing problems with pickups and deliveries. European Journal of Operational Research, 162, 126–141. Ropke, S., & Pisinger, D. (2006). A unified heuristic for a large class of vehicle routing problems with backhauls. European Journal of Operational Research, 171(3), 750–775. Toth, P., & Vigo, D. (2001). The vehicle routing problem, monographs on discrete mathematics and applications. Philadelphia: SIAM. Wolsey, L. A. (1998). Integer programming. New York: John Wiley & Sons. Zachariadis, E. E., Tarantilis, C. D., & Kiranoudis, C. T. (2010). An adaptive memory methodology for the vehicle routing problem with simultaneous pick-ups and deliveries. European Journal of Operational Research, 202, 401–411.