Next Article in Journal
Extending 3D-GIS District Models and BIM-Based Building Models into Computer Gaming Environment for Better Workflow of Cultural Heritage Conservation
Next Article in Special Issue
Finding Effective Item Assignment Plans with Weighted Item Associations Using A Hybrid Genetic Algorithm
Previous Article in Journal
Innovative Hydrogeophysical Approaches as Aids to Assess Hungarian Groundwater Bodies
Previous Article in Special Issue
Performance Evaluation of Hybrid WOA-SVR and HHO-SVR Models with Various Kernels to Predict Factor of Safety for Circular Failure Slope
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

An Attraction Map Framework of a Complex Multi-Echelon Vehicle Routing Problem with Random Walk Analysis

1
Institute of Information Science, University of Miskolc, 3515 Miskolc, Hungary
2
Institute of Logistics, University of Miskolc, 3515 Miskolc, Hungary
*
Author to whom correspondence should be addressed.
Appl. Sci. 2021, 11(5), 2100; https://doi.org/10.3390/app11052100
Submission received: 20 January 2021 / Revised: 19 February 2021 / Accepted: 22 February 2021 / Published: 27 February 2021
(This article belongs to the Collection Heuristic Algorithms in Engineering and Applied Sciences)

Abstract

:
The paper aims to investigate the basin of attraction map of a complex Vehicle Routing Problem with random walk analysis. The Vehicle Routing Problem (VRP) is a common discrete optimization problem in field of logistics. In the case of the base VRP, the positions of one single depot and many customers (which have product demands) are given. The vehicles and their capacity limits are also fixed in the system and the objective function is the minimization of the length of the route. In the literature, many approaches have appeared to simulate the transportation demands. Most of the approaches are using some kind of metaheuristics. Solving the problems with metaheuristics requires exploring the fitness landscape of the optimization problem. The fitness landscape analysis consists of the investigation of the following elements: the set of the possible states, the fitness function and the neighborhood relationship. We use also metaheuristics are used to perform neighborhood discovery depending on the neighborhood interpretation. In this article, the following neighborhood operators are used for the basin of attraction map: 2-opt, Order Crossover (OX), Partially Matched Crossover (PMX), Cycle Crossover (CX). Based on our test results, the 2-opt and Partially Matched Crossover operators are more efficient than the Order Crossover and Cycle Crossovers.

1. Introduction

One of the most important tasks of logistics is the cost-effective delivery of the right goods, to the right place in the right time considering any other constraints. These problem domains belong to the area of the Vehicle Routing Problem (VRP), which has developed in several versions since its first version [1] in 1959. The basic VRP (which is the easiest version of the VRP) is characterized by the fact that the position of a single depot and many customers are known in advance. The demands of the customers are also fixed in advance and, in the system, there is only a single product. The number of vehicles and their capacity constraints is also given in advance. The vehicles start their routes form the depot, they visit some customers, and then return to the depot. The objective function is the minimization of the length of the route. In the following, the main variants of the VRP are detailed that contributed to the development of our own VRP model. Vehicle Routing Problem with Single Depot [2] is a VRP where the system contains a single depot. In case of Vehicle Routing Problem with Multiple Depots [3,4] multiple depots are in the system, and the vehicles start their route from one of the depots, and after visiting the customers, they return to the depot from which they started their route. In the case of the Two-Echelon Vehicle Routing Problem [5], there are not only depots and customers but also satellites. The products are transported from the depot to intermediate locations (satellites) and then from the satellites to the customers. In most of the cases, there are several types of (capacity-constrained) vehicles in the system. If there is only one type of vehicle in the system, it can be said that the system has a homogeneous fleet (Homogeneous Fleet Vehicle Routing Problem [6]). If the system contains several types of vehicles, then the system has a heterogeneous fleet (Heterogeneous Fleet Vehicle Routing Problem [7]). If vehicles have a capacity limit for the goods to be transported, then the system is called Capacitated Vehicle Routing Problem [4,8,9]. Customers to visit can also have a time window, this is called the Vehicle Routing Problem with Time Windows [4].
One of the latest VRP trends is the Fuel-Efficient Green Vehicle Routing Problem [10]; this case prioritizes not the minimization of the length of the route, but the minimization of the fuel consumption of the vehicles as the objective function. In the case of Cumulative Vehicle Routing Problem [11], the minimization of the waiting times for the customers is the objective function. In case of Fuzzy VRP, some factors are given with fuzzy numbers [12]. Milk-run VRP [13] is also a common problem, where the products, materials, etc., are transported between the warehouse and the production line. Another new trend is the application of electric vehicles [14].
Separate algorithms are usually implemented for each VRP model and we also find examples of general heuristics [4]. The investigated system is designed to solve the following VRP model: time windows, capacitated vehicles, multi-depot, and open VRP.
In logistics, beside VRP, we can also find examples of other transportation tasks [15,16], too.
The fitness landscape analysis [17] is a method to investigate the optimization space, it consists of the following elements: the set of possible solution states, neighborhood definition and the fitness function. Although encoding does not belong directly to the concept of the fitness landscape, it plays a big role in the analysis of space. The neighbors of the actual state can be created in several ways (for example, in the case of permutation representation, the 2-opt operator is a common way). To prove that a metaheuristic is suitable for a given optimization task, fitness landscape analysis can be a good alternative. In this case, we compare the main fitness landscape analysis methods with metaheuristics operators. One way of fitness landscape analysis is the trajectory-based sampling approach. In trajectory-based sampling, samples are taken from the search space. We start from a possible state point, then create the state of one or more possible neighbors, then select one of them, which is called the current state. In the literature, the following trajectory-based sampling techniques are investigated [18]: random walk, adaptive walk, reverse adaptive walk, neutral walk, reverse neutral walk, uphill-downhill walk. During a random walk, the next state is the randomly selected neighbor of the current state. In case of an adaptive walk, only a better state is selected from the neighbors of the current state. The reverse adaptive walk is the inverse of the adaptive walk, in which case a state that is worse than the current state is selected from the neighbors of the current state. During a neutral walk, we select a neighbor that shows a fitness value that is close to the current state point. The reverse neutral walk is the opposite of neutral walk, then we select a neighbor where the difference in fitness value is as large as possible. During the uphill downhill walk, we first perform an adaptive step and then a reverse adaptive walk is performed.
This paper focuses on the analysis of the random walk in a complex Vehicle Routing Problem The remainder of the paper is organized as follows: in Section 2, the related works are detailed. In Section 3, the Multi-Echelon Vehicle Routing Problem is presented; in Section 4, the basin of attraction map with random walk is analyzed; in Section 5, the results and discussion is detailed, and after that the summary Section follows.

2. Related Works

In this section, we present some proposals in connection with fitness landscape analysis. Examination of the optimization space is a useful component in analyzing the general behavior of the problem domain. The structure of the optimization space influences the efficiency of convergence, such as the quality of the obtained global or local optimum. Examining the fitness space can be useful for selecting the appropriate optimization methods (which operators are worth using for heuristic algorithms). One benefit of the presented methods is that it can be used as the foundation for further formal theoretical analyzes, there is no need to perform expensive experiments (such as when comparing the results of metaheuristic algorithms with test results).
Our analyzes include Vehicle Routing Problem [19,20], Traveling Salesman Problem [21,22,23], Job Scheduling Problem [21], Test Function [24], Quadratic Assignment Problem [25]. The authors of [19] perform a fitness landscape analysis of the Vehicle Routing Problem. Information is analyzed from a theoretical perspective, and the operators of the Genetic Algorithm (Swap, Inversion, Insertion, Displacement, Partially Matched Crossover—PMX, Uniform Order Crossover—UOX, Cycle Crossover—CX) were chosen as the subject of their analysis. The concepts of information content, partial information content, and density-basin information are reported. The authors of [20] also perform a fitness landscape analysis of the vehicle routing problem. After presenting the mathematical model, the authors review the concepts of search space (fitness landscape). Then, in proportion to the fitness values, the averages of the distances taken from the solutions are illustrated and the distances taken from the best solution are also plotted. The authors of [22] focus on the analysis of the Memetic Algorithm in solving the Traveling Salesman Problem. After discussing each concept of fitness landscape, the authors also plot each solution in a coordinate system according to fitness values. The subject of the analysis is the distance from the optimum and the differences in fitness. Analyses are performed on benchmark data sets. The Traveling Salesman Problem and the No-Wait Job Scheduling are discussed by the authors in [21]. Fitness landscape analysis is performed as autocorrelation, time to local optimum, number of local optimums, its increase as the problem increases, distances from the optimum, the probability that we will reach the optimum. An example of the analysis of benchmark continuous functions can also be found in the literature, for example, in [24] the Differential Evolution Algorithm fitness landscape analysis is presented. Presentations of the most important concepts of fitness landscape: fitness distance correlation, correlation length, epistasis, evolvability and neutrality. The authors focused on a continuous optimization task. Benchmark functions were solved such as Sphere function, Rosenbrock function, Step function, Schwefel function, Rastrigin function, Griewangk function, Ackley function, Easom function, Schwefel’s Ridge function. Schwefel function, the Easom function, the Schwefel’s Ridge function had a small fitness distance correlation, which means, that a Hill Climbing algorithm can easily solve these problems. The analysis of the Traveling Salesman Problem with benchmark data sets is also presented in [23], by comparing the Iterated Local Search and Hill Climbing Algorithm. In this article, experiments were executed on benchmark data sets. They applied 3 TSP-lib instances for the simulations. The distance and running time compared to the best result known so far are compared. An example of the analysis of the Quadratic Assignment Problem (QAP) fitness landscape can also be found in [25]. In addition to presenting the task, operators are also presented as pairwise exchange or swap. In addition, the basin of attraction, the local optima network, the random walk autocorrelation function, and the autocorrelation length are also detailed. The Genetic Algorithm (GA) and Simulated Annealing (SA) were used as algorithms, Partially Matched Crossover (PMX) and Pairwise Exchange Mutation as operators. Correlation analysis is also performed, comparing the two algorithms. The number of local optima and shortest path to the optimum results are also produced. It was concluded that the GA was the stronger algorithm to solve all the studied classes of QAP instances.

3. The Multi-Echelon Vehicle Routing Problem

In this Section our Vehicle Routing Problem model is detailed. Our VRP model and test instance is based on our previous work [26]. Figure 1 illustrates the topology of our problem.
Figure 1 illustrates an example of a VRP sketch on which our own VRP example is based. Figure 1 indicates that our VRP consists of not only customers, indicated with C U = c u 1 , c u 2 , , c u n c and a depot, indicated with D = d 1 , but also satellites. However, the satellites are also different, in levels. The product is shipped from a central depot to the first level satellite S L 1 = s l 1 , 1 , s l 1 , 2 , , s l 1 , n s l 1 , than to the second level satellite S L 2 = s l 2 , 1 , s l 2 , 2 , , s l 2 , n s l 2 , etc., and then to the customers. The system contains a single type of product, indicated with P R = p r 1 . We can also define the set of locations, with the union of the depot, satellites and customers: L O = D S L 1 S L 2 C U = l o 1 , , l o n l o . The system also contains several vehicles, indicated with V E = v e 1 , , v e n v e . Our system contains several components, which can be categorized in the following way: node components, vehicle components, time components, product components, cost components, and metrics (objective function components). To formalize the components, the following notations are introduced: i , j is the node index, k is the vehicle index. The product has not got any special index, because the system contains one type of product. In the following, we use a function-oriented description of the costs associated with the main components of the VRP system in order to provide a more general formalism. The arguments of the functions are equal to the main cost parameters of the VRP.
The objective is the minimization of the length of the route ( min z l r ) , fuel consumption ( min z f c ) , vehicle rental fee ( min z r f ) , route time ( min z r t ) , unvisited customers ( min z u c ) , loading time ( min z l t ) , unloading time ( min z u t ) , loading cost ( min z l c ) , unloading cost ( min z u c ) , administration cost ( min z a c ) , quality control cost ( min z q c c ) and the maximization of the route status ( min z r s ) . For the objective function, it is necessary to normalize the individual components of the objective function. The objective function is a combination of these components, which can be written in the following way: Z = { ( w l r ×   z l r n o r m ) 2 + ( w f c × z f c n o r m ) 2 + w r f × z r f n o r m 2 + w r t × z r t n o r m 2 + w u c × z u c n o r m 2 + w l t   × z l t n o r m 2 + w u t   ×   z u t n o r m 2 + w l c   ×   z l c n o r m 2 + w u c   ×   z u c n o r m 2 + w a c × z a c n o r m 2 + w q c c   ×   z q c c n o r m 2 + w r s   ×   z r s n o r m 2 } 1 / 2 , where z n o r m means the normalization of an objective function component, and w means the weights of each components.
To the node component belongs the travel time between the nodes, which can be written with the following function: T r a v e l T i m e l o i , l o j , v e k , so the travel time varies by locations and vehicle types. The travel distance between the nodes can be written in the following way: T r a v e l D i s t a n c e l o i , l o j , v e k , so the travel distance varies also by locations and vehicle types. Route status between the nodes is also defined: R o u t e S t a t u s l o i , l o j , v e k , it also varies by locations and vehicles. By route condition we mean the technical condition of the route (e.g., potholes) and the security of the route (e.g., robberies). The aim is to maximize the average state value of the generated path. To the vehicle component belongs the capacity constraint of the vehicle, fuel consumption and rental fee. The vehicles have a capacity limit to the product, it can be written in the following way: C a p a c i t y L i m i t v e k , it depends on vehicles and product (but there is only one product in the system). The vehicles have fuel consumption, which varies by vehicles; it can be written in the following way: F u e l C o n s u m p t i o n v e k . In the system, there are own fleet of vehicles and also rented vehicles. The rented vehicles have rental fees, indicated with R e n t a l F e e v e k . The time component indicates components in the system in connection with time, which are the followings: L o a d i n g T i m e l o i , v e k , it indicates, that the loading time depends on locations, and vehicles. In our system there is also unloading time, which can be written in the following way: U n l o a d i n g T i m e l o i , v e k . The loading and unloading time depend on location and vehicle. In the locations, administration time can also arise in connection with the product; the following function indicates this component: A d m i n i s t r a t i o n T i m e l o i , v e k . In connection with product, there is only one component in our system, the product demand of each location, it can be written in the following way: P r o d u c t D e m a n d l o i . The product demand depends on location and product (but there is a single product in our system). Our system consists of the following cost components: loading cost, unloading cost, administrative cost, quality control cost. The loading cost can be written in the following way: L o a d i n g C o s t l o i , v e k and unloading cost can be formulated in the same way: U n l o a d i n g C o s t l o i , v e k . The administration cost is written with A d m i n i s t r a t i o n C o s t l o i , v e k and quality control cost with Q u a l i t y C o n t r o l C o s t l o i , v e k . All of the cost components depend on locations and vehicles.

Our Vehicle Routing Instance

In the tests, we use artificially generated data sets. The test system contains a single depot where the position is generated in interval (0, 100). The satellites in the test system can be divided into two levels. The satellites of the first level were generated with position index in (200, 300), and the positions of the second level are from (400, 500). The number of the satellite nodes is equal to 10 and the number of customers is set to 15. Each level contains two vehicles. The capacity constraints of the vehicles are set to a value between 10,000 and 50,000, and the fuel consumptions are generated between 10 and 100. As detailed above, the costs values in our system are generated from the interval (10, 50). The temporal parameter values are generated in (30, 50). The travel distances between the nodes are calculated with Euclidean distance, the travel time values are from (10, 100), and the route status values are between 100 and 500.

4. Fitness Landscape Analysis: Neighborhood Operators and Random Walk Metrics

Fitness landscape [17] definition covers the set of possible states, the definition of the neighborhood and the fitness function. Many neighborhood operators have appeared in the literature, it depends on whether the problem is a discrete or continuous optimization problem. In the following, the different neighborhood operators of our discrete optimization problem are detailed.
We can distinguish single operand and two-operand operators. In our Multi-Echelon VRP system, the operands are permutations representing the processing order of nodes. During the evaluation, each element of the permutation is lined up and assigned to a vehicle until a solution (vehicle capacity) is encountered. When the constraint is not met, we assign the nodes to a next vehicle.
2-opt [27] operator is widely used in the case of permutation representation of optimization problems (Figure 2).
In the case of 2-opt, two edges are swapped. Swapping means the selection of a given section of the permutation and the reversing of this selected section.
Within the frame of the fitness landscape analysis, we can analyze the performance of the mutation operators. For neighborhood generation the following operators can be used: cycle crossover, order crossover, partially matched crossover.
Figure 3 demonstrates a cycle crossover operator [28]. The cycle crossover locates the cycles in the path. Both parents are permutation, so the cycle crossover takes advantage of the fact that the elements are listed twice (exactly once in both parents). In the example, the (1,9) pair are taken first. The first child collects the first element of the first parent, and the second child receives the first element of the second parent. Then, the second child collects the element denoted by 9. Therefore, we have to look for the item 9 in the first parent, and the first child puts this in the right place. The second child receives the element 1. The first child has already got the element 1, so the circle is closed. The next steps can be summarized in the followings: the first child gets the other elements from the second parent, and the second child gets the other elements from the first parent.
Figure 4 shows the order crossover operator [28]. This operator uses a random fitting section. The position of the chosen fitting section must be the same for both parents. In the second child, the elements of the fitting section of the first parent are substituted with letter H (holes). Similarly, in the first child, the elements of the fitting section of the second parent are substituted with letter H (holes). The holes will be pushed up in the fitting section. After that step, the holes of the first child are substituted with the fitting section of the second parent, and the holes of the second child are substituted with the fitting section of the first parent.
Figure 5 demonstrates the partially matched crossover operator [28]. In the case of this crossover, a random fitting section is chosen. The fitting section of both parents must be in the same position. Pairs are formed from the elements of the fitting section of the two parents. Here, in pairs, the exchange operation itself takes place independently of each other. Pairs denote the elements to be exchanged, but once the elements to be exchanged have been determined, the exchange in both genes and the other gene occurs independently of each other. In our example, we can identify the following pairs: (2,8), (6,5), (7,3), (5,2). We copy parents; these copies will be the new children of the parents. These pairs are changed in the first and then the second child. In our example, we change the 2 with 8, then the 6 with 5, then the 7 with 3 then the 5 with 2.
Within the frame of this operation, it is also important to measure the distance between the potential solutions. In this paper, three techniques are presented, which are the followings: fitness distance, hamming distance, and basic swap sequence distance. Fitness distance is the absolute value of the difference in fitness value between the two solutions. The hamming distance [29] is the number of differences between the representations of the two solutions. The basic swap sequence distance [30] is the number of swaps between the representations of the two solutions.
Basin of attraction are solutions, which convergence to a certain optimum [18]. The goal is to get out of this as easily as possible, because when an algorithm gets to a local optimum, if it has a hard time getting out of it, it has less chance of exploring space, thus reaching the global optimum. Our fitness landscape analysis uses trajectory-based sampling, in which a sample is taken from the search space. In this technique first, a possible state is selected, it will be first the actual state. With the neighborhood operator, we form one or more (in our case one) neighbor states, one of which will be the new current state. In the case of a random walk [18] sampling technique, the random neighbor state of the actual sate will be the actual sate. The neighbor can be better or worse than the actual state. We also prepare an information theory analysis of the fitness landscape, which quantifies the frequency of relative sequences between neighbors. The following important metrics can be evaluated [18]:
  • Partial Information Content: The fitness value of each of the next neighbors is subtracted. It can be a positive number (ascending) that means 1, a negative number (descending), it means −1, or zero (not changing), it means 0. Therefore, the values −1, 0, 1 are indicator values. They show whether the fitness of the neighbor is increasing, decreasing, or remaining of the same. The sequence of the indicator values is given with s o r i g i n a l .
    This procedure removes zero (non-variable) values from the sequence and then leaves only the variable slopes (i.e., descending values after ascending, ascending values after descending), indicated with s n o n z e r o . The value thus obtained is divided by the length of the original sequence.
    M ε = s n o n z e r o s o r i g i n a l
    A value of 1 is rugged while a value of 0 is flat.
  • Expected Partial Information Content: this measure requires the partial information content value. The expected partial information content can be calculated with the following way:
    E M ε = n   ×   M ε 2
    In the formula, n is the length of the s o r i g i n a l , M ε is the partial information content. Expected partial content is not a probability variable, it shows the shape of the fitness landscape.
  • Information Stability: the fitness values of each of the next neighbors are subtracted from each other. The greatest value obtained will be information stability.
  • Entropy: the fitness values of each of the next neighbors are subtracted. It can be a positive number (ascending) that means 1, a negative number (descending), it means −1, or zero (not changing), it means 0. This gives a sequence. Then, the entropy (entropy) is calculated with the following formula, where P i is the relative frequency of the symbols:
    i P i l n P i
  • Information Content: the fitness values of each neighbor solutions are subtracted each other. It can be a positive number (ascending) that means 1, a negative number (descending), it means −1, or zero (not changing), it means 0. This gives a sequence. Then, the information content can be calculated with the following formula, where P p q denotes the relative frequency of the consecutive symbols, where p , q   1 ,   0   , 1 :
    H ε p     q P p q     ×   l o g 6 P p q
  • Density Basin Information: the same as information content, but the basis of the logarithm is different. The Equation of the density basin information is the following:
    h ε p     q P p q   ×   l o g 3 P p q
  • Regularity: the fitness values of each of the next neighbors are subtracted. These values represent regularity.
  • Evolvability portrait: the fitness value of all neighbors of the current solution are divided by the number of neighbors

5. Results and Discussion

In this section, first, the search space is filtered, to find out what solutions space consists of, how they relate to each other, how far apart they are. This step is based on random sampling.
Then, the random walk, which is a trajectory-based sampling, is performed. Results of the random walk are compared with each other, then the result is also analyzed from the point of view of information theory. Operators are good if they give varied results, because then they explore the space well, the basin of attraction is small, so they can easily get out of a local optimum.
The analyzes were performed in Java with our own implementation. The test was running on a personal computer.

5.1. Scattered Filtering of the Search Space

This step means that the search space is sampled. Therefore, we randomly generate solutions and then analyze them. In the tests, we created 100 solutions at random. We chose a relatively low value, even though the task is complicated one, but the number of nodes to be visited is not very large. The following analysis methods were considered: distance of solutions, the distance of the best solution, the distance of the best solution and solutions, cost density. The following distances were analyzed: fitness, hamming, basic swap sequence. For average distances for the search space points were calculated by calculating the average distance to all other search space points. Therefore, each solution is assigned to a point-level average distance.
Figure 6a shows fitness values and average fitness distances. Fitness values range from 120,000 to 160,000 (fitness unit), while average fitness distances range from 5000 to 19,000 (fitness unit). These values also show a parabolic function because small and large fitness values are paired with a large average distance, while medium fitness values are associated with a medium average distance. Figure 6b shows fitness values and average hamming distances. The average hamming distances in each case range from about 35 to 36 (hamming distance unit). Figure 7 shows fitness values and average basic swap sequence distances. The average basic swap sequence distances in each case range from about 28 to 29 (basic swap sequence distance unit). The fitness values and fitness distances from the best solution (the best of the generated solutions) show a linear function, the higher the fitness value of the solution, the greater the distance from the best solution (the best of the generated solutions). According to Figure 6b, the average hamming distances are around 35 (hamming distance unit), while according to Figure 7, the average basic swap sequence distances are around 29 (basic swap sequence distance unit). Fitness values range from 120,000 to 160,000 (fitness unit), and fitness distances range from 4,000 to 36,000 (fitness unit). In case of fitness values and hamming distance from the best solution (the best of the generated solutions), the Hamming distances range from 30 to 40 (hamming distance unit). The basic swap sequence distance from the best solution (the best of the generated solutions) range from 24 to 35 (basic swap sequence distance unit). The fitness distance from the best solution (the best of the generated solutions) and the average of the fitness distances from the other solutions are described by a parabolic function, thus, low and high distances mean high average distances from the best (best of generated solutions) solution, while medium distances mean low average distances from the best (best of generated solutions) solution. The cost density value shows that almost every solution has a different fitness value.
Table 1 presents the meaning of evaluation and Table 2 illustrates the test results.

5.2. Radom Walk Analysis: The Comparison of the Results to Each Other

The average fitness distances (where distances are calculated from the elements of the solution set) during the random walk and 2-opt operator describe a parabolic function based on Figure 8a Fitness values range from 120,000 to 140,000 (fitness unit), while average fitness distances range from 2500 to 9000 (fitness unit). If a solution has a low or high fitness value, the average fitness distance from the other solutions is high, while if it has a medium fitness value, the average distance is low. According to Figure 8b, average hamming distances are not affected by fitness value. Average hamming distances range from 28 to 34 (hamming distance unit). The average basic swap sequence distances (Figure 9) are not affected by fitness value. The average basic swap sequence distances are between 20–26 (basic swap sequence unit).
In the following, the results of order crossover operator is presented, where distances are calculated from the elements of the solution set. The average fitness distances during the random walk and order crossover operators describe a parabolic function. Fitness values range from 120,000 to 140,000 (fitness unit), while average fitness distances range from 2000 to 10,000 (fitness unit). If a solution has a low or high fitness value, the average fitness distance from the other solutions is high, while if it has a medium fitness value, the average distance is low. The average hamming distances are not affected by fitness value. Average hamming distances range from 30 to 34 (hamming distance unit). The average basic swap sequence distances are also not affected by fitness value. The average basic swap sequence distances are between 23–28 (basic swap sequence unit).
In the following, the results of the cycle crossover operator, where the distances are calculated from the elements of the solution set is presented. The distances describe a parabolic function. Fitness values range from 120,000 to 140,000 (fitness unit), while average fitness distances range from 1800 to 5000 (fitness unit). If a solution has a low or high fitness value, the average fitness distance from the other solutions is high, while if it has a medium fitness value, the average distance is low. The average hamming distances are not affected by fitness value. Average hamming distances range from 25 (hamming distance unit) to 30 (hamming distance unit). The average basic swap sequence distances are also not affected by fitness value. The average basic swap sequence distances are between 20–24 (basic swap sequence unit).
In the following, the results of the partially matched crossover operator, where distances are calculated from the elements of the solution set is detailed. The average fitness distances during partially matched crossover operator describes a parabolic function. Fitness values range from 120,000 to 140,000 (fitness unit), while average fitness distances range from 2000 to 7500 (fitness unit). If a solution has a low or high fitness value, the average fitness distance from the other solutions is high, while if it has a medium fitness value, the average distance is low. The average hamming distances are not affected by fitness value. Average hamming distances range from 30 (hamming distance unit) to 34 (hamming distance unit). The average basic swap sequence distances are also not affected by fitness value. The average basic swap sequence distances are between 23–27 (basic swap sequence distance unit).
In the following, the results of the 2-opt operator is presented, where the distances are calculated from the best element of the solution set. Figure 10a shows the fitness values and the fitness distance from the best solution for random walk and 2-opt. The higher the fitness value of a solution, obviously the farther away from the best solution. According to Figure 10b the hamming distances range from 2 to 38 (hamming distance unit). The greater the fitness value of a solution, the greater the hamming distance from the best solution. According to Figure 11, the basic swap sequence distances are between 1–32 (basic swap sequence distance unit). The higher the fitness value of a solution, the greater the distance of the basic swap from the best solution.
In the following, the results of the order crossover operator are detailed, where distances are calculated from the best element in the solution set. The higher the fitness value of a solution, obviously the farther away from the best solution. The hamming distances range from 2 to 40 (hamming distance unit). The greater the fitness value of a solution, the greater the hamming distance from the best solution. The basic swap sequence distances are between 2–34 (basic swap sequence distance unit). The higher the fitness value of a solution, the greater the distance of the basic swap from the best solution.
The results of the cycle crossover operator and a partially matched crossover are detailed, where distances are calculated from the best element of the solution set. These analyzes also show that the higher the fitness value of a solution, the greater the distance of the fitness from the best solution, and the hamming and basic swap sequence distances also increase as a function of the fitness value. For a cycle crossover, the hamming distance is 0–38 (hamming distance unit), while the basic swap sequence distance is 0–29 (basic swap sequence distance unit), and for a partially matched crossover, the hamming distance is 2–40 (hamming distance unit) and the basic swap sequence distance is 1–34 (basic swap sequence distance unit).
In the following, the distance from the best solution and the average of the distances from the other solutions are presented. According to Figure 12a, this is a parabolic function for 2-opt operator, so for large and small distances the average distances are large, while for medium distances the average distances are small for random walk and 2-opt. This shows that the extreme case is not in the middle of the distribution, it is located close to the elements in the center. According to Figure 12b the distance from the best solution during hamming distances does not affect the average distance here. Figure 13 shows the results for the basic swap sequence, here too the average basic swap sequence distance is not affected by the distance from the best solution to the given result. In the case of order crossover, we also get a parabolic function for fitness distances. The results for cycle crossover and partially matched crossover are similar to those for 2-opt and cycle crossover.
Figure 14 and Figure 15 show the cost of density results, i.e., how many solutions with the same fitness value are present between each solution. Based on the results, in the case of 2-opt, almost every solution has a unique fitness value. In the case of order crossover and cycle crossover, some solutions have the same fitness value, while in the case of partially matched crossover, each solution also has a unique fitness value. Table 3. illustrates the evaluation methods of the results of random walk analysis, and Table 4. presents the summary of the results of the random walk where the solutions are compared to each other.

5.3. Random Walk Analysis: Information Content Analysis

In this subsection, the information theory analysis of the random walk technique is evaluated. The following techniques were used: partial information content, expected partial information content, information stability, entropy, information content, density basin information, regularity, and evolvability portrait. Table 5. presents the information content analysis evaluation methods.
The partial information content is 0.55 (partial information content unit) for 2-opt, 0.56 (partial information content unit) for order crossover, 0.47 (partial information content unit) for cycle crossover and 0.63 (partial information content unit) for partially matched crossover. Since the partial information content shows a rugged landscape for 1 and a flat landscape for 0, so in our case we get neither a flat nor a rugged landscape for all neighborhood techniques. The flattest is in the case of a partially matched crossover, so this operator may be the most efficient, while in the case of a cycle crossover it may be the most rugged, so this operator may be the worst.
The expected partial information content values are 27 (expected partial information content unit) for 2-opt, 28 (expected partial information content unit) for OX, 23 (expected partial information content unit) for CX and 31 (expected partial information content unit) for PMX.
The information stability, which shows the largest change in fitness of the neighbors, was given the following values: 7111.15 (information stability unit) for 2-opt, 5809.92 (information stability unit) for OX, 2728.38 (information stability unit) for CX and 5696.12 (information stability unit) for PMX. According to this, the largest change was for 2-opt and the smallest change was for CX.
The lower the entropy, information content, density basin information values, the better, because this means that the search is moving in the same direction (increasing or decreasing the fitness value of the neighbor of the current solutions).
The entropy values are the followings: for 2-opt 1.5809 (entropy unit) for OX 2.5792 (entropy unit) for CX 2.9458 (entropy unit) and for PMX 2.1070 (entropy unit). According to this, the 2-opt has the smallest entropy and the CX has the largest.
The information content values are 0.4101 (information content unit) for 2-opt, 0.8252 (information content unit) for OX, 0.8831 (information content unit) for CX, and 0.6380 (information content unit) for PMX. According to this, the highest value was given by CX and the lowest by 2-opt.
The density basin information values are 0.3285 (density basin information unit) for 2-opt, 0.2814 (density basin information unit) for OX, 0.4183 (density basin information unit) for CX, and 0.2888 (density basin information unit) for PMX.
The diagram of regularity values for each operator (2-opt, OX, CX, PMX) is shown in Figure 16 and Figure 17 which shows ruggedness.
A graph of evolvability portrait values for each operator is shown in Figure 18 and Figure 19, in which we also get a rugged function along the iterations. However, if we plot the values along fitness values, which are shown in Figure 20 and Figure 21. We get a level result (function) that is independent of the fitness value. The higher the value of evolvability, the better the solution, because the more neighbors of the current solution are better than the current solution, so it is more likely to improve during the iteration. In fact, evolvability measures the probability of improvement.

5.4. Analysis of the Test Results

In the case of random walk, the average difference in fitness values of the solutions received by the operators is almost the same, the hamming and basic swap sequence distances are very similar, in the case of order crossover and partially matched crossover the distances are slightly larger than in 2-opt and cycle crossover. Cost of density values also range from 1 to 2, so fitness values are also unique. In the case of a cycle crossover, the cost of density diagram shows that several solutions have the same fitness value, although the number of solutions with the same fitness value is only between 1 and 3. The distances from the best solution also move in large intervals, which also means that the space is well mapped by these operators.
Table 6 analyzes the efficiency of the operators (2-opt, order crossover, cycle crossover, partially matched crossover). The table shows that 2-opt and partially matched crossover (PMX) are more effective than cycle crossover (CX) and order crossover (OX).
The partial information content value was the highest for the partially matched crossover, which means a flat landscape, so this technique proved to be the best during the measurement. In the case of the cycle crossover, the partial information content value is the lowest, so this technique proved to be the worst during the measurement.
The expected partial information content value is also the highest for the partially matched crossover, so this technique proved to be the best.
The information stability value is highest for 2-opt and lowest for cycle crossover. The value shows the largest change in fitness of the neighbors, which means that the higher the value for an operator, the more worthwhile it is to use. In the case of order crossover and partially matched crossover, the value between 2-opt and cycle crossover was obtained, which is almost the same for the two operators.
The lower the entropy, information content, density basin information values, the better the landscape, because this means that the search is moving in the same direction (increasing or decreasing the fitness value of the neighbor of the current solutions).
2-opt has the lowest entropy and CX has the largest. The order crossover and partially matched crossover values are nearly the same, between 2-opt and order crossover.
The 2-opt has the lowest information content value and the cycle crossover has the largest. The order crossover and partially matched crossover values are nearly the same, between 2-opt and order crossover operators.
Of the density basin information values, CX gave the highest value and OX the lowest.
According to the regularity diagram, the difference between the solutions of the cycle crossover is the largest, here we find the largest jumps. In the case of 2-opt, the jumps are the smallest.
Evolvability values are much higher on average for 2-opt and partially matched crossover than for cycle crossover and order crossover, so it is better to use 2-opt and partially matched crossover operators. The higher the value of evolvability, the better the solution, because the more neighbors of the current solution are better than the current solution, so it is more likely to improve during the iteration. In fact, evolvability measures the probability of improvement.

5.5. Verifying our Results with Real Optimization Run

In this subsection, the results of our fitness landscape analysis technique are verified by optimization test runs. The tested operators can also be used in the Genetic Algorithm (GA), so this algorithm was selected for testing. The Genetic Algorithm maintains a population of solutions. Starting from an initial population, it continuously improves elements of the population with crossover and mutation techniques.
Table 7 presents the parameters and results of our test run. The test run was set up with OX in Test 1, PMX in Test 2, CX in Test 3, and 2-opt in Test 4. In Test 5, all three crossing operators were given equal emphasis.
Since the objective is to minimize the fitness value, it can be read from Table 7 that Test 4 was the best, in addition, Tests 1 and 2 also provided high-quality results. The most efficient operator was 2-opt, but the PMX and OX operators also proved to be effective. Using the CX operator and setting equal sharing of each crossover has not been shown to be effective.

6. Conclusions

In this paper, we performed a random walk analysis of the basin of attraction map (which is one of the fitness landscape analyzation methods) for the Multi-Echelon Vehicle Routing Problem. Vehicle Routing Problems is a common optimization task in the field of logistics. The fitness landscape analysis investigates the following components: the set of possible states, the fitness function and the neighborhood. The set of possible states depends on the problem instance and the used operators. The paper presents detailed analysis of the following operators: 2-opt, order crossover, cycle crossover, partially matched crossover. Two investigation methods were used in the random walk analysis. In the first approach, we used the direct solution comparison. Here, we examined the average distances between the solutions, the distances taken from the best solution, and how many solutions have the same fitness value. The following distances were used: fitness distance, hamming distance, basic swap sequence distance. Another approach is the formal information analysis. The following methods were used in the information analysis approach: partial information content, expected partial information content, information stability, entropy, information content, density basin information, regularity, evolvability portrait. Based on the performed tests, the presented fitness landscape analysis method, the basin of attraction map, is suitable for the analysis of the efficiency of some optimization operators. Based on the measurements, the 2-opt and partially matched operators proved to be effective, and the order crossover and cycle crossover proved to be weak.
We can summarize the strengthens of our approach in the followings. It provides an overview of the theoretical foundations and an abstract level analysis of the optimization task. It can be used to explore the search space for the transport task and to select the most efficient operators. The proposed methods can be used also for development and implementation of new heuristic algorithms. Our further research direction is the extension of the analyzes to other heuristic operators.

Author Contributions

Conceptualization, L.K. and A.A.; methodology, L.K., A.A. and T.B.; software, A.A.; validation, L.K. and A.A.; formal analysis, L.K. and A.A.; investigation, L.K. and A.A.; resources, T.B.; data curation, L.K., A.A.; writing—original draft preparation, L.K., A.A. and T.B.; writing—review and editing, L.K., A.A. and T.B.; visualization, L.K., A.A. and T.B.; supervision, L.K. and T.B.; project administration, L.K. and T.B.; funding acquisition, T.B. All authors have read and agreed to the published version of the manuscript.

Funding

The described article was carried out as part of the EFOP-3.6.1-16-2016-00011 “Younger and Renewing University – Innovative Knowledge City—institutional development of the University of Miskolc aiming at intelligent specialisation” project implemented in the framework of the Szechenyi 2020 program. The realization of this project is supported by the European Union, co-financed by the European Social Fund.

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

The data cited in this manuscript are available from the published papers or a corresponding author.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Dantzig, G.B.; Ramser, J.H. The truck dispatching problem. Manag. Sci. 1959, 6, 80–91. [Google Scholar] [CrossRef]
  2. Skrlec, D.; Filipec, M.; Krajcar, S. A heuristic modification of genetic algorithm used for solving the single depot capacited vehicle routing problem. In Proceedings of the Intelligent Information Systems: IIS’97, Grand Bahama Island, Bahamas, 8–10 December 1997; pp. 184–188. [Google Scholar] [CrossRef]
  3. Nagy, G.; Salhi, S. Heuristic algorithms for single and multiple depot vehicle routing problems with pickups and deliveries. Eur. J. Oper. Res. 2005, 162, 126–141. [Google Scholar] [CrossRef] [Green Version]
  4. Pisinger, D.; Ropke, S. A general heuristic for vehicle routing problems. Comput. Oper. Res. 2007, 34, 2403–2435. [Google Scholar] [CrossRef]
  5. Crainic, T.G.; Perboli, G.; Mancini, S.; Tadei, R. Two-echelon vehicle routing problem: A satellite location analysis. Procedia-Soc. Behav. Sci. 2010, 2, 5944–5955. [Google Scholar] [CrossRef] [Green Version]
  6. Ochi, L.S.; Vianna, D.S.; Drummond, L.M.; Victor, A. A parallel evolutionary algorithm for the vehicle routing problem with heterogeneous fleet. Future Gener. Comput. Syst. 1998, 14, 285–292. [Google Scholar] [CrossRef]
  7. Yao, B.; Yu, B.; Hu, P.; Gao, J.; Zhang, M. An improved particle swarm optimization for carton heterogeneous vehicle routing problem with a collection depot. Ann. Oper. Res. 2016, 242, 303–320. [Google Scholar] [CrossRef]
  8. Ralphs, T.K.; Kopman, L.; Pulleyblank, W.R.; Trotter, L.E. On the capacitated vehicle routing problem. Math. Program. 2003, 94, 343–359. [Google Scholar] [CrossRef]
  9. Sitek, P.; Wikarek, J.; Rutczyńska-Wdowiak, K.; Bocewicz, G.; Banaszak, Z. Optimization of capacitated vehicle routing problem with alternative delivery, pick-up and time windows: A modified hybrid approach. Neurocomputing 2021, 423, 670–678. [Google Scholar] [CrossRef]
  10. Poonthalir, G.; Nadarajan, R. A fuel efficient green vehicle routing problem with varying speed constraint (F-GVRP). Expert Syst. Appl. 2018, 100, 131–144. [Google Scholar] [CrossRef]
  11. Kara, İ.; Kara, B.Y.; Yetiş, M.K. Cumulative vehicle routing problems. Veh. Routing Probl. 2008, 1, 85–98. [Google Scholar] [CrossRef] [Green Version]
  12. Bocewicz, G.; Banaszak, Z.; Rudnik, K.; Witczak, M.; Smutnicki, C.; Wikarek, J. Milk-run Routing and Scheduling Subject to Fuzzy Pickup and Delivery Time Constraints: An Ordered Fuzzy Numbers Approach. In Proceedings of the 2020 IEEE International Conference on Fuzzy Systems (FUZZ-IEEE), Glasgow, UK, 19–24 July 2020. [Google Scholar] [CrossRef]
  13. Grzegorz, B.; Izabela, N.; Arkadiusz, G.; Zbigniew, B. Reference model of milk-run traffic systems prototyping. Int. J. Prod. Res. 2020, 1–18. [Google Scholar] [CrossRef]
  14. Zambrano-Martinez, J.L.; Calafate, C.T.; Soler, D.; Lemus-Zúñiga, L.G.; Cano, J.C.; Manzoni, P.; Gayraud, T. Centralized route-management solution for autonomous vehicles in urban areas. Electronics 2019, 8, 722. [Google Scholar] [CrossRef] [Green Version]
  15. Burduk, A.; Musiał, K. Optimization of Chosen Transport Task by Using Generic Algorithms. Lect. Notes Comput. Sci. 2016, 9842, 197–205. [Google Scholar] [CrossRef]
  16. Kłosowski, G.; Gola, A.; Amila, T. Computational Intelligence in Control of AGV Multimodal Systems. Ifac-Pap. 2018, 51, 1421–1427. [Google Scholar] [CrossRef]
  17. Pitzer, E. Applied Fitness Landscape Analysis. Ph.D. Dissertation, Johannes Kepler UniversitÄt Linz, Technisch-Naturwissenschaftliche Fakultät, Linz, Austria, 2013. [Google Scholar]
  18. Pitzer, E.; Affenzeller, M. A comprehensive survey on fitness landscape analysis. In Recent Advances in Intelligent Engineering Systems; Studies in Computational Intelligence; Fodor, J., Klempous, R., Suárez Araujo, C.P., Eds.; Springer: Berlin, Heidelberg, Germany, 2012; p. 378. [Google Scholar] [CrossRef]
  19. Ventresca, M.; Ombuki-Berman, B.; Runka, A. Predicting genetic algorithm performance on the vehicle routing problem using information theoretic landscape measures. Lect. Notes Comput. Sci. 2013, 7832, 214–225. [Google Scholar] [CrossRef]
  20. Czech, Z.J. Statistical measures of a fitness landscape for the vehicle routing problem. In Proceedings of the 2008 IEEE International Symposium on Parallel and Distributed Processing, Miami, FL, USA, 14—18 April 2008; pp. 1–8. [Google Scholar] [CrossRef]
  21. Tayarani-N, M.H.; Prügel-Bennett, A. An analysis of the fitness landscape of travelling salesman problem. Evol. Comput. 2016, 24, 347–384. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  22. Fitness Landscapes and Graphs: Multimodularity, Ruggedness and Neutrality. Available online: http://www-lisic.univ-littoral.fr/~verel/talks/2tut16-verel.pdf (accessed on 12 December 2020).
  23. Fitness Landscapes of Combinatorial Problems and the Performance of Search Algorithms. Available online: https://philippe-preux.github.io/papiers/lil-97-13.pdf (accessed on 12 December 2020).
  24. Uludağ, G.; Uyar, A.Ş. Fitness landscape analysis of differential evolution algorithms. In Proceedings of the 2009 Fifth International Conference on Soft Computing, Computing with Words and Perceptions in System Analysis, Decision and Control, Famagusta, Cyprus, 2–4 September 2009; pp. 1–4. [Google Scholar] [CrossRef]
  25. Chicano, F.; Daolio, F.; Ochoa, G.; Vérel, S.; Tomassini, M.; Alba, E. Local optima networks, landscape autocorrelation and heuristic search performance. Lect. Notes Comput. Sci. 2012, 7492, 337–347. [Google Scholar] [CrossRef] [Green Version]
  26. Kovács, L.; Agárdi, A.; Bányai, T. Fitness Landscape Analysis and Edge Weighting-Based Optimization of Vehicle Routing Problems. Processes 2020, 8, 1363. [Google Scholar] [CrossRef]
  27. Englert, M.; Röglin, H.; Vöcking, B. Worst case and probabilistic analysis of the 2-Opt algorithm for the TSP. Algorithmica 2014, 68, 190–264. [Google Scholar] [CrossRef] [Green Version]
  28. Hussain, A.; Muhammad, Y.S.; Sajid, M.N.; Hussain, I.; Shoukry, A.M.; Gani, S. Genetic algorithm for traveling salesman problem with modified cycle crossover operator. Comput. Intell. Neurosci. 2017, 7430125. [Google Scholar] [CrossRef] [PubMed]
  29. Zhu, K.Q. A diversity-controlling adaptive genetic algorithm for the vehicle routing problem with time windows. In Proceedings of the 15th IEEE International Conference on Tools with Artificial Intelligence, Sacramento, CA, USA, 3–5 November 2003; pp. 176–183. [Google Scholar] [CrossRef]
  30. Wang, K.P.; Huang, L.; Zhou, C.G.; Pang, W. Particle swarm optimization for traveling salesman problem. In Proceedings of the 2003 International Conference on Machine Learning and Cybernetics (IEEE Cat. No.03EX693), Xi’an, China, 5 November 2003; Volume 3, pp. 1583–1585. [Google Scholar] [CrossRef]
Figure 1. Multi-Echelon Vehicle Routing Problem model.
Figure 1. Multi-Echelon Vehicle Routing Problem model.
Applsci 11 02100 g001
Figure 2. 2-opt operator.
Figure 2. 2-opt operator.
Applsci 11 02100 g002
Figure 3. Cycle crossover (CX) operator.
Figure 3. Cycle crossover (CX) operator.
Applsci 11 02100 g003
Figure 4. Order crossover (OX) operator.
Figure 4. Order crossover (OX) operator.
Applsci 11 02100 g004
Figure 5. Partially Matched crossover (PMX) operator.
Figure 5. Partially Matched crossover (PMX) operator.
Applsci 11 02100 g005
Figure 6. The fitness distance (a) and the hamming distance (b) of solutions.
Figure 6. The fitness distance (a) and the hamming distance (b) of solutions.
Applsci 11 02100 g006
Figure 7. The basic swap sequence distance of solutions.
Figure 7. The basic swap sequence distance of solutions.
Applsci 11 02100 g007
Figure 8. The fitness values and their average fitness distances (a) and the average hamming distances (b) from each other in case of 2-opt operator.
Figure 8. The fitness values and their average fitness distances (a) and the average hamming distances (b) from each other in case of 2-opt operator.
Applsci 11 02100 g008
Figure 9. The average basic swap sequence distances from each other in case of 2-opt operator.
Figure 9. The average basic swap sequence distances from each other in case of 2-opt operator.
Applsci 11 02100 g009
Figure 10. The fitness values and their fitness distances (a) and hamming distances (b) from the best of the solutions in case of 2-opt operator.
Figure 10. The fitness values and their fitness distances (a) and hamming distances (b) from the best of the solutions in case of 2-opt operator.
Applsci 11 02100 g010
Figure 11. The basic swap sequence distances from the best of the solutions in case of 2-opt operator.
Figure 11. The basic swap sequence distances from the best of the solutions in case of 2-opt operator.
Applsci 11 02100 g011
Figure 12. The distance from the best solution and the average of the fitness distances (a) and the hamming distances (b) from the other solutions in case of 2-opt operator.
Figure 12. The distance from the best solution and the average of the fitness distances (a) and the hamming distances (b) from the other solutions in case of 2-opt operator.
Applsci 11 02100 g012
Figure 13. The basic swap sequence distances from the other solutions in case of 2-opt operator.
Figure 13. The basic swap sequence distances from the other solutions in case of 2-opt operator.
Applsci 11 02100 g013
Figure 14. The cost of density in case of 2-opt (a) and order crossover operator (b).
Figure 14. The cost of density in case of 2-opt (a) and order crossover operator (b).
Applsci 11 02100 g014
Figure 15. The cost of density in case of cycle crossover (a) and partially matched crossover operator (b).
Figure 15. The cost of density in case of cycle crossover (a) and partially matched crossover operator (b).
Applsci 11 02100 g015
Figure 16. The regularity in case of 2-opt (a) and order crossover operator (b).
Figure 16. The regularity in case of 2-opt (a) and order crossover operator (b).
Applsci 11 02100 g016
Figure 17. The regularity in case of cycle crossover (a) and partially matched crossover operator (b).
Figure 17. The regularity in case of cycle crossover (a) and partially matched crossover operator (b).
Applsci 11 02100 g017
Figure 18. The evolvability along the iterations in case of 2-opt (a) and order crossover operator (b).
Figure 18. The evolvability along the iterations in case of 2-opt (a) and order crossover operator (b).
Applsci 11 02100 g018
Figure 19. The evolvability along the iterations in case of cycle crossover (a) and partially matched crossover operator (b).
Figure 19. The evolvability along the iterations in case of cycle crossover (a) and partially matched crossover operator (b).
Applsci 11 02100 g019
Figure 20. The evolvability in case of 2-opt (a) and order crossover operator (b).
Figure 20. The evolvability in case of 2-opt (a) and order crossover operator (b).
Applsci 11 02100 g020
Figure 21. The evolvability in case of cycle crossover (a) and 2-opt operator (b).
Figure 21. The evolvability in case of cycle crossover (a) and 2-opt operator (b).
Applsci 11 02100 g021
Table 1. Analysis of the scattered filtered search space: meaning of evaluation.
Table 1. Analysis of the scattered filtered search space: meaning of evaluation.
Analysis of Filtered Search Space
TypeOptimal Value
Fitness valuesSmall value difference is optimal. Large value difference between lower and upper bound in case of high mountains and valleys
Average of fitness distances
Average of hamming distancesSmall lower and upper bound are optimal. At high average distances, it is more difficult for the optimization algorithm to navigate through the field because the task is complicated
Average of basic swap sequence distances
Fitness distances of the best solutionA small difference in value is optimal. At great distances, we are far from the optimum
Hamming distances of the best solutionSmall lower and upper bound are optimal. At large distances, it is more difficult for the optimization algorithm to get to the optimum
Basic swap sequence distances of the best solution
Cost densitySmall lower and upper bound are optimal. For high values, the solutions are varied
Table 2. Analysis of the scattered filtered search space: results.
Table 2. Analysis of the scattered filtered search space: results.
Analysis of Filtered Search Space
TypeDistanceLower BoundUpper Bound
Fitness values Fitness (fitness unit)120,000160,000
Average of fitness distances Fitness (fitness unit)500019,000
Average of hamming distances Hamming (hamming distance unit)3536
Average of basic swap sequence distances Basic swap sequence (basic swap sequence distance unit)2829
Fitness distances of the best solution Fitness (fitness unit)400036,000
Hamming distances of the best solutionHamming (hamming distance unit)3040
Basic swap sequence distances of the best solutionBasic swap sequence (basic swap sequence distance unit)2435
Cost density-13
Table 3. The evaluation methods of the results of random walk analysis.
Table 3. The evaluation methods of the results of random walk analysis.
Random Walk Analysis: The Comparison of Results to Each Other
TypeOptimal Value
Fitness valuesIf the value difference between the lower and upper bound is large, the algorithm detects the space well. The lower bound value, on the other hand, must be small
Average of fitness distancesIt explores space well over long distances
Average of hamming distances
Average of basic swap sequence distances
Fitness distances of best solutionLarge upper bound and small lower bound. Then, it explores the space well (large upper bound) but also found a very good solution due to the small upper bound.
Hamming distances of best solution
Basic swap sequence distances of best solution
Cost of densityAt high values, it explores space well
Fitness distance of filtered global optimaLarge upper bound and small lower bound. Then, it explores the space well (large upper bound) but also found a very good solution due to the small upper bound.
Hamming distance of filtered global optima
Basic swap sequence distance of
filtered global optima
Table 4. The summary of the results of the random walk where the solutions are compared to each other.
Table 4. The summary of the results of the random walk where the solutions are compared to each other.
TypeDistanceLower BoundUpper Bound
2-opt
Fitness valuesFitness (fitness unit)120,000140,000
Average of fitness distancesFitness (fitness unit)25009000
Average of hamming distancesHamming (hamming distance unit)2834
Average of basic swap sequence distancesBasic swap sequence (basic swap sequence distance unit)2026
Fitness distances of best solutionFitness (fitness unit)100017,000
Hamming distances of best solutionHamming (hamming distance unit)238
Basic swap sequence distances of best solutionBasic swap sequence (basic swap sequence distance unit)132
Cost of density12
Fitness distance of filtered global optimaFitness (fitness unit)100018,000
Hamming distance of filtered global optimaHamming (hamming distance unit)3040
Basic swap sequence distance of filtered global optimaBasic swap sequence (basic swap sequence distance unit)2234
Order crossover
Fitness valuesFitness (fitness unit)120,000140,000
Average of fitness distancesFitness (fitness unit)200010,000
Average of hamming distancesHamming (hamming distance unit)3034
Average of basic swap sequence distancesBasic swap sequence (basic swap sequence distance unit)2328
Fitness distances of best solutionFitness (fitness unit)100016,000
Hamming distances of best solutionHamming (hamming distance unit)240
Basic swap sequence distances of best solutionBasic swap sequence (basic swap sequence distance unit)234
Cost of density13
Fitness distance of filtered global optimaFitness (fitness unit)011,000
Hamming distance of filtered global optimaHamming (hamming distance unit)2838
Basic swap sequence distance of filtered global optimaBasic swap sequence (basic swap sequence distance unit)2232
Cycle crossover
Fitness valuesFitness (fitness unit)120,000140,000
Average of fitness distancesFitness (fitness unit)18005000
Average of hamming distancesHamming (hamming distance unit)2530
Average of basic swap sequence distancesBasic swap sequence (basic swap sequence distance unit)2024
Fitness distances of best solutionFitness (fitness unit)09500
Hamming distances of best solutionHamming (hamming distance unit)038
Basic swap sequence distances of best solutionBasic swap sequence (basic swap sequence distance unit)029
Cost of density13
Fitness distance of filtered global optimaFitness (fitness unit)150010,500
Hamming distance of filtered global optimaHamming (hamming distance unit)3038
Basic swap sequence distance of filtered global optimaBasic swap sequence (basic swap sequence distance unit)2830
Partially Matched Crossover
Fitness valuesFitness (fitness unit)120,000140,000
Average of fitness distancesFitness (fitness unit)20007500
Average of hamming distancesHamming (hamming distance unit)3034
Average of basic swap sequence distancesBasic swap sequence (basic swap sequence distance unit)2327
Fitness distances of best solutionFitness (fitness unit)013,000
Hamming distances of best solutionHamming (hamming distance unit)240
Basic swap sequence distances of best solutionBasic swap sequence (basic swap sequence distance unit)134
Cost of density12
Fitness distance of filtered global optimaFitness (fitness unit)350017,000
Hamming distance of filtered global optimaHamming (hamming distance unit)3040
Basic swap sequence distance of filtered global optimaBasic swap sequence (basic swap sequence distance unit)2234
Table 5. Information content analysis evaluation methods
Table 5. Information content analysis evaluation methods
Information Content Analysis
MethodOptimal Value
Partial Information ContentThe bigger the better
Expected Partial Information Content
Information Stability
EntropyThe smaller the better
Information Content
Density Basin Information
RegularityHigh values
Evolvability portrait
Table 6. Summary results: basin of attraction map with random walk analysis.
Table 6. Summary results: basin of attraction map with random walk analysis.
Summary Results
Analysis MethodEfficient OperatorWeak Operator
The comparison of results to each other2-opt
OX
PMX
Partial Information ContentPMXCX
Expected Partial Information ContentPMXCX
Information Stability2-optCX
Entropy2-optCX
Information Content2-optCX
Density Basin InformationOXCX
Regularity2-optCX
Evolvability portrait2-opt
PMX
CX
OX
Table 7. Real optimization run with Genetic Algorithm.
Table 7. Real optimization run with Genetic Algorithm.
Test 1Test 2Test 3Test 4Test 5
Run time (sec, average)90.490780.7781.6485.3674.61
Population size6060606060
OX rate90%0%0%10%30%
PMX rate0%90%0%10%30%
CX rate0%0%90%10%30%
2-opt rate (mutation)10%10%10%50%10%
Fitness value (average)119,579119,715124,504117,856128,832
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Agárdi, A.; Kovács, L.; Bányai, T. An Attraction Map Framework of a Complex Multi-Echelon Vehicle Routing Problem with Random Walk Analysis. Appl. Sci. 2021, 11, 2100. https://doi.org/10.3390/app11052100

AMA Style

Agárdi A, Kovács L, Bányai T. An Attraction Map Framework of a Complex Multi-Echelon Vehicle Routing Problem with Random Walk Analysis. Applied Sciences. 2021; 11(5):2100. https://doi.org/10.3390/app11052100

Chicago/Turabian Style

Agárdi, Anita, László Kovács, and Tamás Bányai. 2021. "An Attraction Map Framework of a Complex Multi-Echelon Vehicle Routing Problem with Random Walk Analysis" Applied Sciences 11, no. 5: 2100. https://doi.org/10.3390/app11052100

APA Style

Agárdi, A., Kovács, L., & Bányai, T. (2021). An Attraction Map Framework of a Complex Multi-Echelon Vehicle Routing Problem with Random Walk Analysis. Applied Sciences, 11(5), 2100. https://doi.org/10.3390/app11052100

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop