• No results found

A Hybrid Genetic Algorithm to Solve the Multi-Depot Vehicle Routing Problem

N/A
N/A
Protected

Academic year: 2021

Share "A Hybrid Genetic Algorithm to Solve the Multi-Depot Vehicle Routing Problem"

Copied!
50
0
0

Bezig met laden.... (Bekijk nu de volledige tekst)

Hele tekst

(1)

Multi-Depot Vehicle Routing Problem

Master’s Thesis Operations Research

(2)

Abstract

This thesis investigates a variant of the Multi-Depot Vehicle Routing Problem (MDVRP) in which the sources, supplying the customers, are related. The sources include a distribution center and several intermediate points, where these intermediate points can be physical stores, depots and/or hubs. The relatedness of the sources is due to the assumption that the interme-diate points are supplied by the distribution center, which is one of the sources as well. The objective is to determine a set of routes that minimizes total costs such that capacity constraints and route duration constraints are satisfied and all customers are served. Due to the related-ness of the sources, the cost structure in this variant of the MDVRP is more complex than it is in the regular MDVRP. A hybrid genetic algorithm is presented for the variant of the MDVRP. The algorithm is benchmarked by the nearest neighbour algorithm, which searches for the nearest unvisited customer location until every location is visited, satisfying capacity constraints and route duration constraints. It is shown that the presented algorithm performs well on numer-ous different input sets.

(3)

Contents

1 Introduction 1

1.1 Situation . . . 1

1.2 Motivation . . . 1

1.3 Goal and Scope . . . 2

1.4 Methodology . . . 3

1.5 Structure of Remainder . . . 5

2 Problem Formulation 6 3 Literature Review 11 4 Algorithm 14 4.1 Generation of Initial Population . . . 16

4.1.1 Initial Population Size . . . 16

4.1.2 Assignment Procedure . . . 16 4.1.3 Routing Algorithm . . . 16 4.1.4 Feasibility Check . . . 17 4.2 Evaluation . . . 17 4.3 Education . . . 19 4.4 Generation of Offspring . . . 20

4.4.1 Stopping Criteria Corresponding to Generation of Offspring . . . 20

4.4.2 Parent Selection . . . 21 4.4.3 Create Offspring . . . 21 4.4.4 Routing Algorithm . . . 21 4.4.5 Feasibility Check . . . 21 4.4.6 Evaluation . . . 21 4.4.7 Education . . . 22

4.4.8 Condition 1: Maximum Size Subpopulation . . . 22

4.4.9 Condition 2: Maximum Number of Iterations . . . 22

4.5 Stopping Criterion . . . 22

5 Computational Experiments 23 5.1 Parameter Specification . . . 23

5.2 Experimental Design . . . 24

5.2.1 Relevance of Education Step . . . 25

5.2.2 Influence of Infeasibility . . . 25

5.2.3 Restrictions Due to High-Level of Factor 1 . . . 26

5.2.4 Summarizing the Restrictions . . . 26

(4)

5.2.6 Demand and Service Time Vector . . . 27

5.2.7 Vehicle Capacity and Maximum Route Duration . . . 27

5.2.8 Transport and Handling Costs . . . 27

5.2.9 Stopping Criteria . . . 28

6 Results 29 6.1 Parameter Specification . . . 29

6.1.1 Fluctuating µ and λ . . . 29

6.1.2 Fluctuating nc . . . 34

6.2 Results 2 × 2 × 2 Factorial Experiment . . . 35

(5)

List of Figures

1.1 Direct Distribution . . . 2

1.2 Indirect Distribution . . . 2

2.1 Possible Routing Scheme, Direct Distribution . . . 10

2.2 Possible Routing Scheme, Both Direct and Indirect Distribution . . . 10

4.1 Flowchart General Scheme Algorithm . . . 15

6.1 Input set 1, execution 1 . . . 31

6.2 Input set 1, execution 2 . . . 31

6.3 Input set 2, execution 1 . . . 32

6.4 Input set 2, execution 2 . . . 32

6.5 Input set 3, execution 1 . . . 33

6.6 Input set 3, execution 2 . . . 33

6.7 Parameter set 1: µ = 9, λ = 28 . . . 35

6.8 Parameter set 2: µ = 15, λ = 46 . . . 35

6.9 Parameter set 3: µ = 36, λ = 109 . . . 35

1 Input set 4, execution 1 . . . 44

2 Input set 4, execution 2 . . . 44

3 Input set 5, execution 1 . . . 45

4 Input set 5, execution 2 . . . 45

List of Tables

5.1 Ranges Parameters (Vidal et al., 2011) . . . 23

5.2 Calibration Results (Vidal et al., 2011) . . . 24

5.3 Stopping Conditions (Vidal et al., 2011) . . . 28

6.1 Possible values of nclosefor different nc and µ. . . 34

(6)

Chapter 1

Introduction

1.1

Situation

Nowadays, almost all retailers not only sell their products in physical stores, but also on the internet. This relatively new channel entails several changes for both customers and retail companies. One of these changes concerns the distribution of products to customers. The most classical form of distribution is a product delivery to a physical store, where the customer buys it. However, it is increasingly common that (internet) purchases are delivered at home, especially when there is an additional service included, such as the installation of a product.

1.2

Motivation

As the home delivery of products is increasingly common, it is interesting to investigate the process of distributing products to customers and to examine cost saving possibilities. Retail-ers have the possibility to outsource the distribution (to, in the Netherlands, PostNL or DHL for example) or to execute it themselves. When focusing on retailers that execute the distribution themselves, the most common manner in the Netherlands is direct distribution from a central distribution center (DC) to customers throughout the whole country with vans or trucks. This situation is presented in Figure 1.1. An alternative way to supply customers is through in-direct distribution, where inin-direct implies first distributing products with a larger truck to a store, hub or depot and second, from this intermediate point (IP), supply the home deliveries with smaller vehicles. This situation is presented in Figure 1.2. This indirect distribution might result in cost savings when comparing it to direct distribution. For instance, when considering a retailer with physical stores throughout the country, these stores could be used as interme-diate points for indirect distribution. The stores are provided with usual store supply from the distribution center by large trucks. When using remaining capacity within these trucks to distribute products meant for home deliveries to the stores, i.e. the intermediate points, the products might be closer to the end customers without making (hardly any) extra costs. Notice that a combination of direct and indirect distribution is possible as well: certain home deliveries are supplied directly from a distribution center while others go through intermediate points.

(7)

Figure 1.1: Direct Distribution Figure 1.2: Indirect Distribution

In the literature, the problem of determining how to distribute products when there are mul-tiple possible sources (distribution centers and intermediate points; stores, depots, hubs) is known as the Multi-Depot Vehicle Routing Problem (MDVRP). This is an extension of the clas-sical Vehicle Routing Problem (VRP). The clasclas-sical VRP involves the determination of a set of routes that minimizes total cost or travel distance, while satisfying certain constraints and serving all customers from a single source. The MDVRP minimizes total cost or travel distance while considering multiple sources serving the customers.

The difference between the regular MDVRP and our situation is that the sources are supposed to be unrelated in the regular MDVRP while in our situation they are not. The intermediate points in our situation of the MDVRP are supplied with stock from a distribution center which is one of the sources as well. Using the remaining capacity of the vehicles to supply the interme-diate points with products meant for home deliveries, results in a more complex cost structure and hence, the objective function of minimizing costs is more complex.

The regular MDVRP is a well known problem in the literature and nowadays it is still ex-amined a lot. This literature will be discussed later on and the application of the literature to cover our situation of the MDVRP will be presented as well.

1.3

Goal and Scope

The goal of this research is to solve the variant of the MDVRP with related sources. We want to determine a near optimal routing scheme for the daily distribution of all products with home delivery service when assuming that every home delivery can either come from a central dis-tribution center or from an intermediate store, hub or depot. This goal can be summarized as follows:

Determine a model that specifies how to distribute home delivery orders, when multiple sources supplying these orders are available, while minimizing the total cost.

(8)

Let us sum up the following assumptions to create a clear view on the situation under consid-eration, i.e. to define the scope of this research:

• Home delivery from an intermediate point is always possible, i.e. the storage capacity and inventory level at the intermediate points are always large enough.

• The capacity of the trucks used to supply the intermediate points is unrestricted, i.e. there is always enough remaining capacity in the trucks to distribute home delivery orders to the intermediate points.

• Two types of costs corresponding to the distribution of the products are included: trans-port costs, which are presented as a cost per kilometer, and handling costs at the interme-diate points, which are presented as a fixed cost per unit of volume. There are no costs included corresponding to the transport of products from the distribution center to an intermediate point, since existing transport lanes are used, i.e. as it is assumed that the capacity of the trucks corresponding to these existing transport lanes is unrestricted, no additional costs are included. Hence, the cost per kilometer corresponding to the trans-port only concerns the actual routes starting and ending at the sources.

• Orders are handled on a daily basis. Consolidation of orders over time is not considered. • Return flows and return handling costs are outside the scope of this research.

• Delivery and departure time windows are outside the scope of this research.

• Vehicles start and end their route at the same source (a distribution center or intermediate point).

• The amount of vehicles available at the sources is unrestricted. When the maximum vehicle capacity is reached, a new route is started and a new vehicle for the new route is assumed to be available.

As the scope of the research is defined, let us now discuss the approach to tackle the research question.

1.4

Methodology

To solve the special situation of the MDVPR, we present a hybrid genetic algorithm. A hybrid algorithm combines several heuristics to exploit their best properties. The genetic algorithm concept, introduced by Holland (1975), is based on mimicking an evolutionary process. The procedure starts with an initial population of candidate solutions where each solution in the population is called a chromosome. In the MDVRP a chromosome represents a routing scheme. Over time the population evolves as iteration by iteration the chromosomes evolve. These iter-ations are called generiter-ations. During each iteration, i.e. during each generation, every chromo-some in the population is evaluated by its fitness measure. This fitness measure represents the quality of the chromosome. In the MDVRP it represents the quality of the routing scheme. This quality might be expressed in the costs related to the routing scheme, but as we will discuss later on, several other aspects may be considered as measurements as well. Chromosomes are selected to execute genetic operations as crossover and mutation according to their fitness. The fitter the chromosome, the higher the probability of being selected which enhances evolution.

(9)

In the crossover phase two chromosomes, called parents, are combined to generate offspring, i.e. to generate a new routing scheme in the MDVRP. The mutation phase is concerned with mutating the offspring to maintain a diverse population. Maintaining a diverse population helps to avoid being trapped in a local optimum which is one of the nice properties of genetic algorithms. When a prespecified number of offspring is created, a new generation is obtained by selecting the best parents and offspring according to their fitness and removing the other chromosomes. After a prespecified number of generations is executed the algorithm returns the best chromosome. This chromosome, hopefully, represents a near optimal solution, i.e. a near optimal routing scheme.

The genetic algorithm we propose to solve the MDVRP determines the routing and assignment scheme for inserted order data. The assignment scheme presents the assignment of the home deliveries to the sources and the routing scheme presents the actual routes. As the number of vehicles available at the sources is assumed to be unrestricted, the hybrid genetic algorithm also presents the required number of vehicles at each source.

The proposed algorithm is based on the heuristic introduced by Vidal et al. (2011). This so called Hybrid Genetic Search with Adaptive Diversity Control (HGSADC) solves three exten-sions of the VRP; the MDVRP, the periodic VRP (PVRP) and the multi-depot periodic VRP (MD-PVRP). Vidal et al. (2011) demonstrate that “for all currently available benchmark instances for the three problem classes, HGSADC identifies either the best know solutions, including the op-timal ones, or new best solutions”. There is no recent research known to improve the method proposed by Vidal et al. (2011). Therefore, this is the method we use to develop our algorithm. The main difference between the heuristic proposed by Vidal et al. (2011) and our heuristic is associated with the relatedness of the sources in our situation due to the distribution center, being a source itself, supplying the other sources. To deal with this relatedness, the handling costs at the intermediate points are included in our algorithm. This implies that home deliv-eries supplied directly from a distribution center are only concerned with a cost related to the actual transportation, while the home deliveries supplied from an intermediate point are not only concerned with a cost related to the transportation but also with a cost related to the extra handling required at the intermediate point.

Another difference between the two heuristics concerns the routing algorithm that is used. Vidal et al. (2011) apply the giant-tour representation and the Split algorithm as Prins (2004) did. We apply a routing algorithm that is easier to build within the used software (MATLAB), namely the nearest neighbour algorithm. As this algorithm is less advanced, the outputted routes are less advanced as well. This problem is overcome by the route improvement phase that is included in the so called education phase of the algorithm.

(10)

The algorithm we propose is tested by a 2 × 2 × 2 factorial experiment. Hence, eight different experiments are investigated to develop a clear overview of the capabilities of the algorithm. The three factors correspond to the total number of locations included (source locations plus customer locations), the demand of the customers and the required service time at the customer locations.

1.5

Structure of Remainder

The next chapter discusses the problem formulation with a formal statement of the MDVRP. Chapter 3 exhibits the Literature Review in which an overview is presented of heuristics de-veloped to solve the MDVRP, including the heuristic of Vidal et al. (2011). In Chapter 4 our hybrid genetic algorithm is described. Chapter 5 discusses the set up of the 2 × 2 × 2 factorial experiment in more detail. The results corresponding to the testing of the algorithm are ana-lyzed in the sixth chapter and conclusions upon these results are drawn in Chapter 7. Finally, in Chapter 8 recommendations for further research are discussed.

(11)

Chapter 2

Problem Formulation

As previously mentioned, the situation we consider can be described by the Multi-Depot Vehi-cle Routing Problem (MDVRP) which is an extension of the classical VehiVehi-cle Routing Problem (VRP). The VRP is a well known problem in several fields of study of operations research due to its application in many industries, such as garbage collection, goods distribution and mail delivery (Liong et al., 2008). Dantzig and Ramser (1959) defined the VRP as a generalization of the Traveling Salesman Problem (TSP). The TSP is concerned with the task of determining the shortest route or least cost route to visit a given set of customer locations exactly once, starting and ending at the same depot. Dantzig and Ramser (1959) extended the TSP by specifying de-mand at each customer location and imposing a capacity constraint on the vehicles supplying these locations. Later on numerous extensions of the VRP have been investigated, including the MDVRP.

Before presenting the formal statement of the MDVRP, let us present an overview of the vari-ables used in this section and the rest of this thesis:

n The total number of customers.

s The total number of sources.

VC The set representing the customers.

VD The set representing the sources.

V = VC∪ VD The set representing all locations together.

qj j ∈ VC The demand of customer j, in a prespecified unit of volume (m2, m3

or pallets).

τj j ∈ VC The service time required by customer j, in minutes.

tij i, j ∈ V The time required to travel from i to j, in minutes.

dij i, j ∈ V The travel distance from i to j, in kilometers.

Q The vehicle capacity, in a prespecified unit of volume (m2, m3or pallets).

T The maximum route duration, in minutes.

c1 The transport costs per kilometer, in euros.

c2 The handling costs per prespecified unit of volume (m2, m3or pallets),

in euros.

R The set representing the separate routes of one routingscheme.

xijk i, j ∈ V, k ∈ R The binary variable describing whether arc (i, j) is included in route k.

zij i, j ∈ V The binary variable describing whether customer j is assigned to

(12)

Now, let us present the mathematical representation of the MDVRP. Let G = (V, A) represent a complete graph where the set V represents the vertices and A the arcs. Let V = VD∪ VC where

VD corresponds to the set containing the sources and VC to the set representing the customers.

Define the number of sources to be s and the number of customers to be n, i.e. VD = {1, ..., s},

VC = {s + 1, ..., s + n} and hence, V = {1, ..., s, s + 1, ..., s + n}. Let the demand and service

duration of customer j ∈ VC be defined by qj and τj respectively. Let Q define the capacity

of the vehicles that distribute the products to the customers and let T represent the maximum amount of time a route can last. The arcs (i, j) ∈ A indicate the travel possibility from vertex i to j, i, j ∈ V . Let the travel distance from i to j, i, j ∈ V be defined by dij. The time required to

travel from i to j is represented by tij, i, j ∈ V . The total duration of a route is the time required

for traveling plus the required service time. Solving the MDVRP involves the determination of a set of routes that minimizes the total duration while satisfying the route duration and vehicle capacity constraint and serving all customers. Another well known objective is minimizing total costs instead of time.

In our special situation of the MDVRP we deal with the objective of minimizing total costs. As discussed in the previous chapter, we include two types of costs: transport costs (a cost per kilometer) and handling costs at intermediate points (a cost per unit of volume). Let the cost per kilometer and the cost per unit of volume be defined by c1and c2respectively. To be able to

define the mathematical representation of our situation of the MDVRP the following variables are required as well as the previously defined variables:

• Let the set R represent the separate routes of one routing scheme.

• Let binary variable xijk, i, j ∈ V and k ∈ R, describe whether arc (i, j) is included in route

kin the routing scheme; hence, xijkis equal to one if vertex j is entered after vertex i in

route k and zero otherwise.

• Let binary variable zij, i ∈ VD and j ∈ VC, describe whether customer j is assigned to

source i; i.e. zij is equal to one if customer j is assigned to source i and zero otherwise.

Now, the costs can be described as follows: transport costs =X i∈V X j∈V X k∈R xijk· dij· c1 (2.1) handling costs = X i∈VD X j∈VC zij · qj· c2 (2.2) total costs = (2.1) + (2.2)

(13)

s.t. X i∈V X j∈V xijk· (tij+ τj) ≤ T ∀k ∈ R (2.3) X i∈V X j∈VC xijk· qj ≤ Q ∀k ∈ R (2.4) X i∈VD X k∈R xijk = 1 ∀j ∈ VC (2.5)

xijk ∈ {0, 1} ∀k ∈ R and ∀i, j ∈ V (2.6)

zij ∈ {0, 1} ∀i, j ∈ V (2.7)

Constraint (2.3) corresponds to the route duration constraint, i.e. every route k ∈ R can last at most T minutes. Hence, ∀k ∈ R we must have τ (k) ≤ T where τ (k) represents the total route duration of route k. This total route duration is the sum of the required travel time and the service duration, i.e. τ (k) =P

i∈V

P

j∈V xijk· (tij + τj). The second constraint represents

the vehicle capacity constraint, i.e. every route k ∈ R can have a maximum load equal to Q. Constraint (2.5) indicates that every customer j ∈ VC is visited exactly once. The fourth and

fifth constraint define the binary variables xijkand zij respectively.

To be able to give more insight into our situation of the MDVRP, consider the following fic-titious example (i.e. no real data is used):

• Assume the existence of one distribution center (DC), defined by 1;

• Assume the existence of two physical stores, which can be used as intermediate points for indirect distribution. Let these intermediate points be defined by 2 and 3;

• Assume the existence of four customers, represented by 4, 5, 6, 7;

From this it follows that we have VD = {1, 2, 3}, VC = {4, 5, 6, 7}and hence, V = VD ∪ VC =

{1, ..., 7}.

Define the following variables for the situation with three sources and four customers, where again no real data is used, i.e. it is a fictitious example:

• Let the following matrix define the travel distances from i to j, i, j ∈ V , in kilometers:

dij =           0 3 5 6 2 5 7.5 3 0 9 3 4.5 7 11 5 9 0 11.5 6 7 3 6 3 11.5 0 8 6.5 12.5 2 4.5 6 8 0 8 9 5 7 7 6.5 8 0 6 7.5 10 3 12.5 9 6 0           ;

(14)

• Let the following vector define the demand of customer j ∈ VC, in m2:

qj =



3 4 2.5 3  ;

• Let the following vector define the required service time of customer j ∈ VC, in minutes:

τj =



30 40 25 30  ;

• Let the maximum amount of time a route can last be equal to T = 180 minutes; • Let the capacity per vehicle be equal to Q = 8 m2;

• Let the transport costs be c1 = 2euro per kilometer;

• Let the handling costs at the intermediate points be c2= 1.50euro per m2.

As previously mentioned, the most common way of supplying customers is through direct dis-tribution, i.e. assigning every customer to the DC. A possible routing scheme corresponding to this assignment scheme is presented in Figure 2.1. This routing scheme consists of two routes: 1) 1 - 4 - 5 - 1 and 2) 1 - 6 - 7 - 1. The first route is 16 kilometers long and lasts 166 minutes. The second route is 18.5 kilometers long and also lasts 166 minutes. The vehicle load corresponding to the routes is 7 m2 and 5.5 m2 respectively. The costs corresponding to this routing scheme

are equal to 2 · (16 + 18.5) = 69 euros. Note, no handling costs are included since there are no deliveries from intermediate points.

Now impose the possibility of indirect distribution. Notice that supplying customer 4 from intermediate point 2 and customer 7 from intermediate point 3 would be a reasonable option to consider. A possible routing scheme corresponding to this situation is presented in Figure 2.2. This routing scheme consists of three routes: 1) 1 5 6 1, 2) 2 4 2 and 3) 3 -7 - 3. The first route length is 15 kilometers and the second and third route are both 6 kilo-meters long. The route durations are 155, 66 and 66 minutes respectively. The vehicle loads are 6.6, 3 and 3 m2 respectively. The total costs corresponding to this routing scheme are 2 · (15 + 6 + 6) + 1.50 · (3 + 3) = 63euros, where 54 euros of these costs are related to the actual transport costs and 9 euros to the handling of the products at the intermediate points. Comparing the two situations indicates that the possibility of indirect distribution results in cost savings.

(15)

Figure 2.1: Possible Routing Scheme, Direct Distribution

Figure 2.2: Possible Routing Scheme, Both Direct and Indirect Distribution

Notice that the routing schemes presented above are two examples of numerous possibilities to solve the considered situation. This situation is quite abstract and far from reality as no real data is used and hence, it is not a surprise that cost savings occur. Although the situation is unrealistic, it illustrates our problem well.

(16)

Chapter 3

Literature Review

Numerous papers have been written on the VRP and several extensions of the problem have been discussed over the years, including the MDVRP. Other extensions include the VRP with time windows, the VRP with backhauls, the VRP with pick-ups and deliveries and the VRP with multiple use of vehicles (Gendreau et al., 2008). Nowadays still a lot of research is per-formed on extensions of the VRP. These extensions become increasingly broad. Aforemen-tioned extensions of the VRP are combined and various constraints are added.

Exact algorithms as well as heuristics have been developed to solve the VRP and its numerous extensions. The exact algorithms are time consuming and hardly ever applicable to problems concerned with more than 50 customers (Liong et al., 2008). Therefore when considering real-istic situations heurreal-istics are usually more applicable.

An extremely broad review would be required to discuss the entire range of exact algorithms and heuristics corresponding to the classical VRP and its extensions. Since we are dealing with the MDVRP, the remainder of this literature review focuses on the heuristics developed to solve this problem.

The first heuristics to solve the MDVRP were developed, among others, by Tillman and Cain (1972) and Wren and Holliday (1972). Tillman and Cain (1972) introduced a heuristic based on the savings algorithm proposed by Clarke and Wright (1964) which implies that combining routes might result in cost savings. This heuristic, as many other heuristics developed in the early literature, does not include an improvement stage on the initial solution to obtain a bet-ter solution (Zhang et al., 2011). The heuristic introduced by Wren and Holliday (1972) does include an improvement stage; in the first stage an initial feasible solution is obtained and in the second stage several refinement heuristics are applied to this solution. This procedure of obtaining an initial solution first and then improving it, is one of the two common structures of the heuristics developed to solve the MDVRP (Ho et al., 2008).

(17)

the distances to the depots closest and second closest to the customers. Raft (1982) developed an algorithm for which the problem is decomposed in five, instead of two, stages. The five stages are solved separately and then connected iteratively.

Chao et al. (1993) applied the cluster first - route second principle presented by Golden et al. (1977) to obtain an initial solution and improve this solution by changing the depot assignments of the customers. Hence, they combine the two commonly applied structures in heuristics solv-ing the MDVRP. In the improvement procedure Chao et al. (1993) applied the record-to-record approach (Dueck, 1993). To explain this record-to-record approach define the best solution ob-tained so far as S, an alternative solution as S0, the total distance corresponding to solution S as record R and deviation D as a preset percentage of R. The alternative solution S0is selected as new solution if the objective value of S0 is less than R + D. Hence, the record-to-record ap-proach allows deteriorations to prevent from being trapped in a local optimum to eventually obtain a good solution.

Tabu search algorithms on the MDVRP have been first introduced by Renaud et al. (1996) and Cordeau et al. (1997). These algorithms use memory structures to memorize the examined so-lutions. A tabu search algorithm will not examine a solution more than once as it marks the visited solutions as tabu in its memory for a certain number of iterations (Glover, 1990). During the 21st century research has primarily been dedicated to extensions of the MDVRP but the classical MDVRP has been investigated as well. Pisinger and Ropke (2007) developed a unified heuristic to solve five extensions of the classical VRP including the MDVRP. This heuris-tic uses the adaptive large neighbourhood search framework designed by (Ropke and Pisinger, 2006). This framework expands and contracts the search for a better solution by choosing adap-tively among a number of heuristics that are able to insert and remove customers in a route. Ho et al. (2008) were the first to propose two hybrid genetic algorithms to solve the MDVRP. The first algorithm generates the initial population randomly while the second algorithm uses the savings algorithm and the nearest neighbour heuristic. The nearest neighbour heuristic is a routing heuristic that solves the TSP. It searches for the nearest unvisited customer lo-cation until every lolo-cation is visited and then returns to the starting lolo-cation. Mirabi et al. (2010) presented three hybrid heuristics to solve the MDVRP where these heuristics combine components from constructive heuristic search and improvement methods. With constructive heuristic search a complete solution, i.e. all customers satisfied, is ’constructed’ by extending an empty solution. The three hybrid heuristics presented by Mirabi et al. (2010) differ in the improvement technique they apply.

(18)

schemes” of including feasible as well as infeasible solutions and refreshing the population when a prespecified maximum number of iterations made without improving the best solution obtained so far is reached. Also, a diversity contribution is included in the fitness measure used to evaluate the chromosomes. This enhances the evolution of the population while avoiding it to converge prematurely.

Another difference is the inclusion of an education phase to further enhance the improvement of chromosomes. In this phase local-search procedures educate the chromosomes by improv-ing their routimprov-ing schemes.

As previously mentioned, Vidal et al. (2011) demonstrated that their algorithm is the best per-forming method to approach the MDVRP. Therefore, this is the method we use to develop our own algorithm which is presented in the next chapter.

(19)

Chapter 4

Algorithm

The general scheme of our proposed meta-heuristic is presented below in text, General Scheme Algorithm, as well as in a flowchart, Figure 4.1. The heuristic generates an initial population of feasible and infeasible chromosomes which are separated in subpopulations. These chro-mosomes represent a routing scheme corresponding to a randomly determined assignment scheme. Next, the chromosomes undergo some route improvement techniques during the so called education step. Thereafter, the roulette wheel selection operation (Golberg, 1989) is ap-plied to select parents which assignment schemes are randomly combined to produce offspring and a routing scheme corresponding to this new created assignment scheme is determined. The created offspring undergoes the education procedure and is included in the subpopulation ac-cording to its feasibility.

General Scheme Algorithm

Step 0: While the computation time & the amount of iterations allowed without improve-ment both have not reached their maximum do:

Step 1: Generate initial population.

1.1 Whilethe maximum size of the initial population is not reached do:

1.2Obtain an assignment scheme by randomly assigning every customer to a source.

1.3Determine a routing scheme corresponding to the assignment scheme by applying the near-est neighbour algorithm.

1.4 Ifthe route duration constraint is not satisfied by one of the routes from the routing scheme obtained in 1.3 then insert the chromosome into the infeasible subpopulation. Else insert the chromosome into the feasible subpopulation.

1.5 End while.

Step 2: Evaluate every chromosome by its evaluation value.

Obtain the evaluation value for every chromosome in the population where this is a combina-tion of the rank of the costs of the chromosome and the rank of the diversity contribucombina-tion of the chromosome.

Step 3: Educate the chromosomes.

The chromosomes are educated by applying some techniques to improve their routing schemes.

Step 4: Generate offspring, i.e. new chromosomes.

(20)

number of iterations allowed in a row without improving the best solution obtained so far is not reached do:

4.2Select parents according to the roulette wheel selection operation (Golberg, 1989).

4.3Apply random crossover to the assignment schemes of the parents to produce offspring.

4.4Determine a routing scheme corresponding to the assignment scheme of the offspring by applying the nearest neighbour algorithm.

4.5 Ifthe route duration constraint is not satisfied by one of the routes from the routing scheme obtained in 4.4 then insert the offspring into the infeasible subpopulation. Else insert the off-spring into the feasible subpopulation.

4.6Obtain the evaluation value of the offspring.

4.7Educate the offspring.

4.8 Ifone of the two subpopulations reached its maximum size then select survivors and return to 4.1.

4.9 Ifthe prespecified number of iterations allowed in a row without improving the best solu-tion obtained so far is reached then diversify the populasolu-tion by selecting survivors and return to 1.1.

4.10 End while. Step 5: End while.

Figure 4.1: Flowchart General Scheme Algorithm

(21)

In the following sections the algorithm will be explained in more detail. The numbering of the sections and subsections is similar as the numbering of the steps within the General Scheme Algorithm. The first section considers the generation of the initial population. Second, the evaluation process of the chromosomes is described. The third section is concerned with the process corresponding to the education of the chromosomes. Next, the generation of offspring is discussed. Last, the fifth section discusses the stopping criterion. As previously mentioned, our heuristic is based on the method presented by Vidal et al. (2011). We will indicate the differences between this heuristic and ours. When no difference is indicated, one can assume similarity between the two methods.

4.1

Generation of Initial Population

4.1.1 Initial Population Size

In the first step of the algorithm, an initial population is generated containing 4µ chromosomes, where µ is the minimum size of both the feasible and infeasible subpopulation. The set P op represents the chromosomes. At the end of the initialization step, at least one of the subpop-ulations contains µ chromosomes. The other one might contain less than µ chromosomes and hence, is then incomplete since µ is defined as the minimum size of both subpopulations. The maximum size of each subpopulation is µ + λ. A subpopulation reaching its maximum size undergoes a procedure of selecting µ survivors. This survivor selection procedure is discussed in subsection 4.4.8.

Let us, in the next subsections explain the actual generation of the initial population.

4.1.2 Assignment Procedure

Chromosomes are created by randomly assigning the n customer locations to a source i where i represents the distribution center, stores and/or depots, i = 1, ..., s. Hence, we are dealing with ncustomers locations and s sources; therefore we can define the following sets of locations VD = {1, ..., s}, VC = {s + 1, ..., s + n}and V = VD∪ VC = {1, ..., s, s + 1, ..., s + n}. The

assign-ment of customer location j ∈ VC, to source i ∈ VD, can be defined by yj = iand hence, the

as-signment scheme of chromosome r is represented by an array of n asas-signments: [ys+1, ..., ys+n].

This assignment scheme for chromosome r can also be defined by matrix Z(r) := zij(r)where

zij(r) = 1if customer location j ∈ VC is assigned to source i ∈ VD, and zero otherwise. Both

yj = iand Z(r) are decision variables representing the same decision. 4.1.3 Routing Algorithm

The actual chromosomes not only represent an assignment scheme, but also the routes obtained by applying the nearest neighbour heuristic to this assignment scheme. Here, our algorithm differs from the one presented by Vidal et al. (2011) as they apply the giant-tour representation and the Split algorithm as Prins (2004) did.

As previously described, the nearest neighbour heuristic is a routing heuristic that solves the TSP. It searches for the nearest unvisited customer location until every location is visited and then returns to the starting location. In our case the starting location is always one of the sources i ∈ VD. We apply an extension of the classical nearest neighbour heuristic by inserting

a capacity constraint. From source i ∈ VDthe nearest neighbour j ∈ VC, for which yj = iholds,

(22)

entered. This continues until the maximum vehicle capacity is reached or all customer locations j ∈ VC, for which yj = iholds, are visited. When the maximum vehicle capacity is reached

while there remain unvisited customers, the route is ended by returning to i and a new route is started. After applying the nearest neighbour heuristic, the matrix X(r) := xijk(r)represents

the routing scheme corresponding to chromosome r. Hence, it describes which arcs (i, j) are driven in chromosome r in route k ∈ R(r), i.e. xijk(r) = 1if arc (i, j) is driven in k ∈ R(r) and

zero otherwise where the set R(r) represent the routes corresponding to chromosome r.

4.1.4 Feasibility Check

Next, the feasibility of the chromosomes is investigated. When the route duration constraint is not satisfied by one of the routes in the routing scheme of the chromosome, it is inserted in the infeasible subpopulation. When every route satisfies the route duration constraint, the chromosome is inserted in the feasible subpopulation. In mathematics this yields when τ (k) > T for at least one k ∈ R(r), where τ (k) represents the route duration of route k and T the maximum route duration, the solution is infeasible and hence, it is inserted in the infeasible subpopulation; otherwise it is inserted in the feasible subpopulation.

4.2

Evaluation

The chromosomes are evaluated by their fitness measure. From now on this fitness measure is represented by the so called evaluation function. This function combines the rank of the costs of the chromosome and the rank of the diversity contribution of the chromosome. As Vidal et al. (2011) describe, this construction of evaluating the individuals not only by their objective value, i.e. their costs, which is most common, but by also including a diversity contribution show to “not only efficiently avoid premature population convergence, but also outperforms traditional diversity management methods relative to the general behavior of the solution method”. The costs are composed of two components: total costs (i.e. transport and handling costs) and a penalty cost corresponding to the feasibility of the chromosome. The two ranks correspond to the location of the chromosome in the ordered list of costs and diversity contributions.

Let us discuss the parts of the evaluation function separately in both words as well as in math-ematics in the following subsections.

Total Costs

For the total costs corresponding to the routing scheme of the chromosome r we distinguish two costs: the transport costs and the handling costs. The transport costs are the costs related to the actual driving from A to B. We define c1to be a fixed cost per kilometer. In this cost per

kilometer multiple factors influencing the transport costs are taken into account such as fuel costs, driver costs, insurances, fleet costs. Now, in mathematics we can define the transport costs for chromosome r to be:

ctransport(r) = X i∈V X j∈V xijk(r) · dij · c1,

where dij represents the distance of driving from i to j. These distances are acquired using the

program XCargo.

Define the handling costs by c2, a cost per unit of volume (per m2, per m3, per pallet). The

(23)

quantity used is dependent on the investigated data. Let demand at customer location j ∈ VC

be defined by qj, where qj is defined in the same quantity. In mathetics we have:

chandling(r) = X i∈VD X j∈VC zij(r) · qj· c2.

The total costs for chromosome r are represented by the sum of the transport and handling costs, i.e. ctotal(r) = ctransport(r) + chandling(r).

Penalized Costs

The penalty for infeasibility is the amount of time that exceeds the route duration constraint multiplied by a penalty parameter, ω. This parameter dynamically changes during the execu-tion of the algorithm to increase the favorability of naturally feasible soluexecu-tions. The penalty pa-rameter is initially set equal to 1 and is adjusted according to the number of iterations executed without improving the best solution obtained so far. This adjustment can be mathematically defined by:

ω = 1 + ITactual ITallowed

,

where ITactual is the actual number of iterations executed without improvement so far and

ITallowedthe maximum allowed number of iterations without improving the best solution

ob-tained so far. Notice, ω ∈ [1, 2]. Hence, during the execution of the algorithm, the penalty parameter increases and hence, the penalty costs weigh more.

In mathematics, the penalized cost for chromosome r is defined by: φ(r) = ω · X

k∈R(r)

max {0, τ (k) − T } .

Diversity Contribution

The diversity contribution of chromosome r is represented by ∆(r). It is “the average distance to its ncloseclosest neighbours” (Vidal et al., 2011). The average distance is a figurative distance;

it represents the average difference between the assignment schemes of the individuals. Let the following equation represent the mathematic notation of the average distance:

γ(r1, r2) = 1 n X j∈VC [I(δj(r1) 6= δj(r2))]

where I(δj(r1) 6= δj(r2)) is an indicator function returning one if the source assignment of

customer location j ∈ N is not the same for chromosomes r1and r2 and zero otherwise.

The diversity contribution is now mathematically represented by: ∆(r) = 1

nclose

X

r2∈Nclose

γ(r, r2)

(24)

Evaluation Function

Now let the actual evaluation function be defined by: EV (r) = ξc· RAN K  1 ctotal(r) + φ(r)  + ξdc· RAN K(∆(r)). RAN K  1 ctotal(r)+φ(r) 

and RAN K(∆(r)) present the rank of individual r corresponding to the costs and the diversity contribution respectively. The rank of the costs considers the inverse of the sum of the total costs and penalized costs to be able to assign the highest rank to the chro-mosome with the lowest costs. In this situation the best chrochro-mosome r maximizes the evalua-tion value EV (r). The ranks are weighted by the parameters ξcand ξdc. These parameters are

initially equal to 0.5 and are dynamically adjusted during the execution of the algorithm. The adjustments can be mathematically presented in the following way:

ξc= 0.5 + 0.5 ITactual ITallowed ξdc= 0.5 − 0.5 ITactual ITallowed

The adjustments make sure that during the execution, the costs become heavier weighted while the weight on the diversity contribution becomes lighter.

4.3

Education

When the evaluation process is complete, the 4µ initial chromosomes are educated. This ed-ucation procedure is concerned with nine route improvement techniques. To discuss these techniques first define the following variables:

• Let u be a customer location, i.e. u ∈ VC;

• Let route k(u) represent the route containing customer location u in chromosome r, k(u) ∈ R(r);

• Let (u, v) be the partial route from u to v in k(u);

• Let the neighbourhood of u be defined by the h · n, h ∈ [0, 1], closest neighbours of u; • Let v be a customer location in the neighbourhood of u;

• Let x be the successor of u in k(u) and let y be the successor of v in k(v).

For the nine route improvement techniques the following restrictions are required: u 6= y and v 6= x. Now the following nine moves can be defined:

1. Remove u and place it after v;

2. Remove (u, x) and place (u, x) after v; 3. Remove (u, x) and place (x, u) after v; 4. Swap u and v;

(25)

5. Swap (u, x) and v; 6. Swap (u, x) and (v, y);

7. If k(u) = k(v), replace (u, x) and (v, y) by (u, v) and (x, y); 8. If k(u) 6= k(v), replace (u, x) and (v, y) by (u, v) and (x, y); 9. If k(u) 6= k(v), replace (u, x) and (v, y) by (u, y) and (x, v).

The last two moves are inter-route moves, the seventh move is an intra-route move, the first three moves are so called insertions and moves four to six are swaps.

During the education step every chromosome is educated which, as the route improvement techniques suggest, yields improving the routing schemes obtained by the nearest neighbour algorithm. A customer location u is randomly selected as well as customer location v from the neighbourhood of u. The nine moves are randomly investigated and the first move yielding an improvement of the evaluation value, EV (r), is implemented. Then, for the same chromo-some, a new u and v are selected and again the nine moves are randomly investigated. This process of randomly selecting u and v continues until a u is selected for which none of the nine moves results in improvement. When this event occurs, the feasibility of the chromosome is in-vestigated and the new evaluation value of the chromosome, EV (r), is obtained. This implies that the feasibility check is performed numerous times during the education step as well as the evaluation phase. When the feasibility of the chromosome changed, it is deleted from the old subpopulation and inserted in the other one. Thereafter, the process starts all over again for the next chromosome until all chromosomes have been educated.

4.4

Generation of Offspring

4.4.1 Stopping Criteria Corresponding to Generation of Offspring

The initial population is extended by creating offspring while two conditions are satisfied: • the two subpopulations have not reached their maximum size, (condition 1),

• and the prespecified number of iterations allowed in a row without improving the best solution obtained so far is not reached (condition 2).

When one of the two subpopulations reaches its maximum size a prespecified number of sur-vivors is selected. With this group of sursur-vivors the procedure of generating offspring starts all over again.This process is explained in more detail in section 4.4.8.

Let the number of iterations allowed in a row without improving the best solution obtained so far be defined by L = 0.75 · λ. When the number of iterations in a row without improving the best solution obtained so far reaches L, the population is diversified by first selecting a prespecified number of survivors and then returning to the first step of the algorithm to create a new initial population containing the selected survivors. This diversification process is dis-cussed in more detail in subsection 4.4.9

(26)

4.4.2 Parent Selection

The first step in creating offspring is the selection of parents. We will apply the roulette wheel selection operation (Golberg, 1989) which takes the evaluation of the chromosomes into account to select the best chromosome with the highest probability and the worst with the lowest. Here our algorithm differs from the one presented by Vidal et al. (2011) as they randomly select par-ents.

The roulette wheel selection operation, as the name suggests, is based on a roulette wheel. Let the interval [0, 1] represent the roulette wheel. Every chromosome is represented on the roulette wheel, i.e. on the interval [0, 1], by a specific part. Let the probability of obtaining chromosome r be defined by:

P (r) = P EV (r)

∀r∈P opEV (r)

.

Now the part on the interval [0, 1] related to chromosome r can be defined by :   X i∈{1,...,r−1} P (i), X i∈{1,...,r} P (i)  .

To illustrate this consider the following:

• chromosome r = 1 corresponds to the part [0, P (1)],

• chromosome r = 2 corresponds to the part [P (1), P (1) + P (2)],

• chromosome r = 3 corresponds to the part [P (1) + P (2), P (1) + P (2) + P (3)], • etcetera.

To select two parents two random numbers on the interval [0, 1] are selected; i.e. the roulette wheel is spinned twice. The chromosomes corresponding to these two randomly generated numbers are selected as parents.

4.4.3 Create Offspring

To create offspring the assignment schemes of the two selected parents are randomly combined. This random crossover results in a new assignment scheme for a new chromosome. From this point onwards the newly created offspring follows the same procedure as the chromosomes do in the process of generating the initial population:

4.4.4 Routing Algorithm

• the nearest neighbour algorithm is applied to obtain a routing scheme corresponding to the assignment scheme of the new chromosome,

4.4.5 Feasibility Check

• its feasibility is investigated,

4.4.6 Evaluation

(27)

4.4.7 Education

• the chromosome is educated and inserted in the subpopulation corresponding to its fea-sibility.

4.4.8 Condition 1: Maximum Size Subpopulation

As previously described, offspring is created as long as the following two conditions are satis-fied:

• the two subpopulations have not reached their maximum size µ + λ,

• and the prespecified number of iterations allowed in a row without improving the best solution obtained so far, L, is not reached.

In this section we will describe what happens when one of the two subpopulations reaches its maximum size. In the next section we will discuss what happens when the second condition is not satisfied.

When one of the two subpopulations reaches its maximum size, µ + λ, the µ best chromosomes are selected for both subpopulations. These two groups of µ chromosomes are the survivors. The two remaining groups, both containing λ chromosomes, are deleted from the subpopula-tions. With the two groups of µ chromosomes the whole process re-enters the stage of generat-ing offsprgenerat-ing until one of the two conditions is unsatisfied again.

4.4.9 Condition 2: Maximum Number of Iterations

When the number of iterations in a row without improving the best solution obtained so far reaches L, the µ/3 best chromosomes are selected for both subpopulations. These two groups of µ/3 chromosomes are the survivors. The remaining two groups of chromosomes are deleted from the subpopulations. With the two groups of µ/3 chromosomes the whole process re-enters the stage of creating an initial population to diversify the chromosomes. The population is expanded until it has size 4 · µ again and the process continues with this new initial population.

4.5

Stopping Criterion

The whole process evolves generation by generation until either one of the two stopping crite-ria is met. The two stopping critecrite-ria are a maximum computation time, Tmax, and a maximum

number of iterations the algorithm is allowed to make without improving the best solution ob-tained so far, ITallowed. When the computation time reaches Tmax or the number of iterations

ITactualreaches ITallowed, the algorithm is stopped and the best chromosome, i.e. the best

rout-ing scheme, obtained so far is returned as the best solution.

(28)

Chapter 5

Computational Experiments

The hybrid genetic algorithm presented in the previous chapter is tested by a 2 × 2 × 2 factorial experiment. The factors correspond to the total number of locations included (source locations plus customer locations), the demand of the customers and the required service time at the customer locations. For each factor two levels are specified: a high and a low level. Before discussing these levels in more detail, the parameter specification is investigated.

5.1

Parameter Specification

In the previous chapter the following four parameters have been defined: µ, λ, nclose and h.

Vidal et al. (2011) specified the following ranges for these parameters to use in the calibration of their values:

Parameter Range

µ Minimum population size [5,200]

λ Number of offspring in a generation [1,200]

nc Proportion of close individuals considered in the diversity contribution, [0,0.25] nclose = nc · µ

h Proportion of close individuals considered in the education step [0,1] Table 5.1: Ranges Parameters (Vidal et al., 2011)

According to Vidal et al. (2011) these ranges are appropriate due to values found in the litera-ture (this holds for the subpopulation sizes), conceptual requirements (“a local distance mea-sure is assumed to implicate not more than 25% of the population”) and parameter definition (this holds for the last proportion). Due to programming reasons we need to deal with the fol-lowing two constraints: λ must be larger than 3 · µ and µ must be divisible by three. Therefore we define the ranges [6, 66] and [19, 199] for µ and λ respectively. For nc and h we apply the ranges Vidal et al. (2011) set.

(29)

Parameter Value

µ Minimum population size 25

λ Number of offspring in a generation 70

nc Proportion of close individuals considered in the diversity contribution, 0.2 nclose= nc · µ

h Proportion of close individuals considered in the education step 0.4 Table 5.2: Calibration Results (Vidal et al., 2011)

Instead of applying the CMA-ES as Vidal et al. (2011) did, we apply a trial and error procedure to obtain parameter values within the previously mentioned ranges. The results according to the specification of the parameters will be presented in the next chapter.

When the parameter values are specified, the algorithm can be tested. Let us, in the next sec-tion, discuss the experimental design we apply.

5.2

Experimental Design

As previously mentioned, a 2 × 2 × 2 factorial experiment is applied. Before discussing the factors in more detail, let us present an overview of the input the algorithm requires:

1. A distance matrix, including the distances dij, i, j ∈ V , representing the distance from i

to j in kilometers;

2. A time matrix is required, including the travel times tij, i, j ∈ V , representing the time

required to travel from i to j in minutes;

3. A demand vector, including the demands qj, j ∈ VC, representing the demand of

cus-tomer j in m2, m3or in pallets (depending on the investigated situation);

4. A service time vector, including the service time τj, j ∈ VC, representing the service time

required at customer j in minutes;

5. The vehicle capacity, Q (in m2, m3 or in pallets; again depending on the investigated situation);

6. The maximum amount of time (in minutes) a route can last, T ; 7. The transport costs c1, per kilometer;

8. The handling costs at the intermediate points, c2.

Some of these requirements are defined by the factors of the experiment. Let us therefore now discuss these factors. For each factor two initial levels (high and low) are specified. We mention initial since preliminary experiments demonstrated that the possible levels are restricted. Be-fore discussing the results corresponding to the preliminary experiments and the consequences for the possible levels, let us discuss the levels we initially defined:

Factor 1:The number of locations included (source locations plus customer locations).

(30)

• Low: 50 zip codes including one distribution center, 4 intermediate points and 40 cus-tomer locations.

Factor 2:The demand of the customers.

A customer’s demand is randomly assigned from a set of three values. We assume that orders are not piled up during transport and therefore define demand in m2.

• High: {4; 5, 5; 5.5} m2;

• Low: {0.5; 1; 1.5} m2.

Factor 3:The required service time at the customer locations.

A customer’s service time is randomly assigned from a set of three service times defined in minutes.

• High: {60; 90; 120} minutes; • Low: {10; 30; 60} minutes.

As previously mentioned, preliminary experiments indicated that the algorithm encounters some issues when dealing with these initial levels. Let us therefore summarize the findings corresponding to the preliminary experiments in the following subsections and indicate the consequences for the algorithm and the levels of the factors.

5.2.1 Relevance of Education Step

One of the findings does not include any problems and does not have consequences for the levels of the factors, but weighs heavily on the execution of the algorithm and is therefore mentioned first. The education step accounts for a large part of the total computation time and does not influence the eventual best solution, independent of the parameter settings and independent of the input. Therefore, it is excluded from the algorithm which implies that h does not need to be defined.

5.2.2 Influence of Infeasibility

A problem was encountered corresponding to the infeasibility of solutions. A solution is infea-sible when the route duration constraint is exceeded by at least one route within the solution. Independent of the parameter settings and independent of the input, the appearance of infea-sible solutions is enlarged by the assignment procedure we apply. As customers are randomly assigned to a source instead of to a source nearby, the distances between the locations within one route are most of the time quite large, especially when the set of locations is small. Hence, the time required to travel is quite large as well, which results in a large route duration. The two extreme settings of high- and low-levels of the factors, both increase the appearance of infeasibility. In the low-levels setting this is due to the fact that the demands within the de-mand range are small and in the high-levels setting this is due to the fact that the service times within the service time range are large. In the low-levels setting this can be explained as fol-lows: the routing algorithm, i.e. the nearest neighbour algorithm, inserts customers into a route as long as the vehicle capacity allows this. This implies that when demand of all customers is relatively small, many customers are included in a single route which results in a large route duration and hence, many infeasible solutions. Large service times, in the high-levels setting,

(31)

in combination with large traveling times due to the random assignment, also results in large route durations.

The actual problem we encountered is the fact that for the two extreme settings only infea-sible solutions were included. Starting with only infeainfea-sible solutions, the chances of including feasible solutions during the execution of the algorithm are rare. Therefore, to be able to in-vestigate the extreme situations the decision was made to change the assignment procedure to increase the appearance of feasible solutions. Instead of only assigning customers randomly, from now on customers are with 70% chance assigned randomly and with 30% chance they are assigned to one of their three nearest sources. Assigning customers to one of their nearest sources, decreases the traveling times and therefore increases the appearence of feasiblity. This does not solve the entire problem. The low-level demand range {0, 5; 1; 1, 5} m2still allows

too many customers to be in one route. Even when all customers are assigned to their nearest source, the vehicle capacity is so large, that many customers are included in one route and the route duration constraint is exceeded. Therefore, from now on the low-level demand range is specified as {1.5; 2; 2.5} m2. Then, independent of the parameter settings, not only infeasible

solutions are included when investigating the low-levels setting.

For the high-levels setting the problem is also not entirely solved by the new assignment pro-cedure. Therefore, the new high-level service time range is specified as {30; 60; 90} minutes. Then, independent of the parameter settings, feasible solutions are included as well.

5.2.3 Restrictions Due to High-Level of Factor 1

For the high-level setting, we initially intended to investigate 500 zip codes. Unfortunately MATLAB encounters memory issues when dealing with such a large number of locations. Therefore, the high-level of Factor 1 is adjusted to 350 zip codes.

MATLAB is able to cope with 350 zip codes, although it still encounters memory issues when dealing with large values of µ and λ, i.e. large populations. The maximum values for µ and λ MATLAB is able to deal with are 15 and 46 respectively. This implies that for the high-levels setting new ranges for µ and λ need to be considered, namely [6, 15] and [19, 46] respectively. For the low-level setting the larger, previously defined ranges are investigated.

5.2.4 Summarizing the Restrictions

Results according to the preliminary experiments have changed the levels of the factors within the 2 × 2 × 2 factorial experiment quite drastically. Let us therefore present a clear overview of the new levels to be able to discuss the results corresponding to the different experiments in the next chapter.

Factor 1:The number of locations included (source locations plus customer locations).

• High: 350 zip codes including one distribution center, 9 intermediate points and 340 cus-tomer locations;

(32)

Factor 2:The demand of the customers.

• High: {4.5; 5; 5.5} m2;

• Low: {1.5; 2; 2.5} m2.

Factor 3:The required service time at the customer locations.

• High: {30; 60; 90} minutes; • Low: {10; 30; 60} minutes.

Combining these three factors, with each two levels, results in eight different experiments. These experiments are all investigated to develop an overview of the capabilities of the algo-rithm. The results corresponding to the investigation of the experiments will be presented in the following chapter. However, before discussing the results, let us return to the input require-ments which will be discussed in the following sections.

5.2.5 Distance and Time Matrix

The first two input requirements, a distance and time matrix, are dependent on Factor 1: The number of locations included. The by Factor 1 defined number of zip codes are randomly drawn from the total of 4039 possible four-digit zip codes from The Netherlands. Next, the program XCargo is used to obtain both the distance matrix as well as the time matrix.

5.2.6 Demand and Service Time Vector

The demand vector and service time vector are dependent on the level of Factor 2 and Factor 3 respectively. As previously mentioned, the demand vector is obtained by randomly assigning demand from the by Factor 2 specified set to each customer. The service time vector is obtained by randomly assigning service time from the by Factor 2 specified level to each customer.

5.2.7 Vehicle Capacity and Maximum Route Duration

As previously mentioned, due to the assumption that orders are not piled up during transport, we use the measurement m2 to define the vehicle capacity. Expert knowledge is used to define

the capacity to be 15 m2and a maximum route duration of 10 hours, i.e. 600 minutes. 5.2.8 Transport and Handling Costs

The program NEA is used to obtain a transport cost per kilometer ofe1.53. In this cost several components are included such as fuel costs, driver costs, insurances, fleet costs and mainte-nance costs.

Expert knowledge is used to define the average time required to load and unload a vehicle and the cost per hour for an employee to do this work. This information is combined to acquire the handling costs ofe4.25 per m2.

(33)

5.2.9 Stopping Criteria

Besides the input requirements, the stopping criteria need to be defined as well.

Vidal et al. (2011) applied three different sets of stopping conditions for ITallowed and Tmax,

namely:

(ITallowed, Tmax)

(104, 10min) (2 · 104, 30min) (5 · 104, 60min)

Table 5.3: Stopping Conditions (Vidal et al., 2011)

We initially intend to investigate these same sets of stopping criteria. Whether the algorithm is fast enough to indeed apply these sets, will be discussed in the following chapter.

(34)

Chapter 6

Results

This chapter contains two sections of results according to two sets of experiments. The first section is concerned with the results according to the parameter specification and the second deals with the results obtained when executing the algorithm with the, in the first section, specified parameter values. Throughout this whole chapter the defined stopping criterion is ITallowed= 104iterations.

6.1

Parameter Specification

To identify appropriate parameter values, several different input sets are investigated. For each input set, the extreme parameter values, i.e. the minima and maxima within their ranges, are investigated first. Next, a trial and error procedure is applied to obtain the most suitable values for µ and λ. Let us summarize our findings in the following subsections. Notice that the parameter h is not investigated, since the education step is excluded from the algorithm.

6.1.1 Fluctuating µ and λ

Unfortunately it is hard to draw conclusions upon the most suitable setting of the parameter values µ and λ. Varying the settings of the parameter values for numerous input sets does not contribute to unambiguous results. This is illustrated by discussing several different input sets, but before doing so, let us present three statements that hold in general, independent of the parameter settings and input sets:

• The smaller µ and λ, the faster the algorithm reaches the stopping criterion, i.e. the smaller the required execution time.

• The larger µ and λ, the smaller the number of generations included while executing the algorithm.

(35)

Although the drop down occurs faster for larger µ and λ, the drop down may not be as large as a drop down occurring for a smaller µ and λ, as illustrated by the figures. Hence, it is not possible to state that a larger µ and λ are more suitable.

Let us investigate different input sets to illustrate that we are not capable of drawing strict con-clusions upon the most suitable settings of the parameters µ and λ. Each input set is investigate for a fixed value of nc, namely nc = 0.20. The following section discusses the actual influence of fluctuating the parameter value nc.

(36)

Input set 1

The first investigated input set includes 25 locations: one distribution center, two intermediate points and 22 customer locations. The demand range is [2; 3; 4] m2and the service time range is

[50; 60; 70]minutes. Figures 6.1 and 6.2 present results for several examined settings for µ and λ. For only 25 locations included, small µ and λ are sufficient, i.e. include better solutions, and therefore larger settings for µ and λ are excluded from the figures. The two figures correspond to obtained solutions of executing the algorithm twice for the exact same input set. Notice that the patterns in the two figures are quite similar. For the different settings of µ and λ the algorithm converges to approximately the same solution. As the algorithm is fastest for the smallest µ and λ, and the obtained final solution does not differ for different settings of the parameter values, the figures suggest that the setting with the smallest value is best to apply, i.e. µ = 6 and λ = 19. Notice that these two figures illustrate only two examples of several executions of the algorithm for different demand and service time ranges. The observed patterns hold in general. Hence, for an input set including 25 customers a µ and λ equal to 6 and 19 respectively, are suitable, independent of the demand and service time range.

Figure 6.1: Input set 1, execution 1

Figure 6.2: Input set 1, execution 2

(37)

Input set 2

Figures 6.3 and 6.4 illustrate results obtained by executing the algorithm twice for the same input set. This input set includes 50 locations: one distribution center, four intermediate points and 45 customer locations. The demand range is [2; 3.5; 5] m2 and the service time range is [45; 60; 75]minutes. Again, large settings for µ and λ are excluded since small parameter values are most suitable. Notice that the results corresponding to the two execution are quite similar for µ = 6 and λ = 19, although the patterns corresponding to the other settings of µ and λ differ quite a lot. Other executions of the algorithm for the same input set yield again other patterns. This implies that it is not possible to specify the most suitable values of µ and λ for this specific input set. Executing the algorithm for other input sets containing 50 locations, also results in numerous observed patterns and no unambiguous specification of the most suitable parameter values.

Figure 6.3: Input set 2, execution 1

(38)

Input set 3

The third investigated input set includes 200 locations: one distribution center, nine intermedi-ate points and 190 customer locations. The demand range is [3; 3.5; 4] m2and the service time

range is [10; 20; 30] minutes. Figures 6.5 and 6.6 present the results corresponding to several settings of parameter values. Notice that the overall patterns within the two figures are quite similar, although looking closer suggests that it is again hard to determine the most suitable parameter values for µ and λ. More executions of the algorithm for both this input set as well as for other input sets containing 200 locations, imply this same result of being incapable of setting µ and λ.

Figure 6.5: Input set 3, execution 1

Figure 6.6: Input set 3, execution 2

Referenties

GERELATEERDE DOCUMENTEN

Second, the problem is difficult to solve because of the high degree of freedom in its variables (Agra et al. The problem of routing and scheduling is related to the vehicle

The vehicle routing problem which is the subject of this paper is to determine a set of tours of minimum total length such that each of a set of customers, represented by points in

If in this situation the maximum number of reduced rest periods are already taken, while a split rest of 3 hours together with the customer service time still fits within the 15

Since we show that the DDVRP is a generalization of the Multi Depot Vehicle Routing Problem (MDVRP), we can benchmark the ALNS solutions against best known solutions to

In this paper we study the two echelon vehicle routing problem with covering options (2E-VRP- CO), which arises when both satellite locations and covering locations are available

For larger input instances the GA was again benchmarked against both specified practical distri- bution rules. These experiments showed that the performance of the GA decreases in

This is true since it is the best solution in all solution spaces cut away by the piercing cuts and the remaining part of the solution space cannot contain a better solution, since

We propose a multi-start simheuristic that combines biased-randomized versions of classical routing and packing heuristics for solving the introduced stochastic variant of the