# Multidepot Two-Echelon Vehicle Routing Problem for Earthwork Allocation Optimization

Authors: Qinglong Zhang; Naifu Deng; Yanwen Zhu; Zhenping Huang
Journal: Mathematical Problems in Engineering
 (2022)

---

## Abstract

Prior to the construction of most engineering projects, earthwork is a complex and time-consuming task, requiring iterative operations in civil engineering. The effectiveness of earthworks determines the cost of many AEC (architecture, engineering, and construction) projects (e.g., road, embankment, railway, and slope engineering). As a result, creating effective earthwork planning is critical. The earthwork allocation problem is simplified in this study to the vehicle route problem (VRP), which is often studied in the field of transportation and logistics. An optimization model for the earthwork allocation path based on the modified genetic algorithm with a self-adaptive mechanism is developed to work out the global optimal hauling path for earthwork. The findings of the study are also used to shape the basic topographic shape of the Winter Olympic Skiing Course Project. Furthermore, a comparative study with the former methods is conducted to validate the performance of our proposed method on tackling such a multidepot two-echelon vehicle routing problem. Because of its flexibility, this optimization model is extremely compatible with various evolutionary methods in many fields, making future development viable and practicable.


---

## Body

## 1. Introduction

Earthwork allocation is a substantial and repetitive task required for the majority of AEC (architecture, engineering, and construction) projects [1, 2]. Earthwork is the process of leveling or shaping the ground in a target area by moving or handling the geological materials that make up the target area. This geological material allocation usually includes excavation, loading, handling, unloading, and compaction operations in different areas. Sometimes it may also include some intermediate steps, such as material mixing or processing [3, 4]. Because of the complexity of earthwork activities, the volume of work is large and often accounts for more than 50% of the total project cost [5]. In particular, among the above steps, transporting geotechnical materials from one location to another is in most cases the most expensive [6, 7]. Since earth distribution is usually a continuous and repeatable iterative process, finding an effective technical approach to rationalize the earth deployment path can yield significant economic benefits, such as reduced fuel consumption or carbon emissions [1, 8].Earthwork allocation is essentially an extended application of the vehicle path problem (VRP). VRP is a combinatorial dynamic planning problem that seeks to utilize a fleet of vehicles, such as trucks, to serve a certain number of customers with different cargo requirements under certain constraints, such as the delivery time or load capacity of nominated vehicles [9, 10]. It primarily refers to a set of issues in which an ideal route for a fleet of vehicles based on one or more depots may be identified for a number of geographically separated consumers [11, 12]. Figure 1 represents the schematic diagram of a single-depot VRP. The vehicle route problem, first proposed by Dantzig and Ramser in 1959 [11], is a critical problem in transportation and logistics. Since its inception, the problem has piqued the interest of experts and researchers, and it has aided the development of transportation, aviation, navigation, communications, electric power, computer science, and other fields in some ways. The VRP, on the other hand, is a nondeterministic polynomial hard (NP-hard) problem that cannot be solved in polynomial time. In other words, given the existing knowledge, it is impossible to find the exact answer to the problem, or it will take ages to find the exact solution. It is especially true when the task is large-scale. As a result, the current research has changed from finding the perfect answer using an exact algorithm to finding a better solution using a heuristic approach.Figure 1
Schematic diagram of a single-depot VRP (also see in [6]).The frequently used heuristic algorithms are the genetic algorithm (GA) [13–16], ant colony algorithm (ACA) [17–19], particle swarm optimization algorithm (PSOA) [20–23], simulated annealing algorithm (SAA) [24–26], tabu search algorithm (TSA) [27–29], differential evolution algorithm (DEA) [30], and so on. ACA is adept at path planning problems, in particular, ACA shows great robustness in solving the traveling salesman problem [31]. However, ACA requires a large amount of computation because it usually requires all ants to choose the same route, which is the optimal line. In practical calculation, it is difficult to achieve this situation under a given number of cycles [32]. SAA is universal and robust, which is suitable for parallel processing and complex nonlinear optimization problem, nevertheless, it also relies on higher computation resources and longer computing time [33]. PSOA and TSA are the opposite, it has fast convergence speed, however, it is easy to produce premature convergence problem, especially in dealing with complex multipeak search problems [33, 34]. DEA and GA are evolutionary algorithms. DEA is modified from GA by improving the mutation operation to accelerate its convergence speed [35], however, fast convergence speed may lead to a premature problem when the initial population is small [36].In this paper, GA is adopted to study the earthwork allocation path problem as it is one of the earliest algorithms to be applied to the field of transportation [37]. Furthermore, research demonstrates that GA is well-suited to handling NP-hard issues, including various variables, parameters, objectives, and weak connectivity across different areas [38, 39]. GA can also record numerous solutions at once, and the simultaneous optimization process for multiple solutions can be used to solve multiobjective optimization problems [39]. Furthermore, GA has a better universality and compatibility, which allows for further improvement when combined with other heuristic algorithms [40].Therefore, this paper focuses on the vehicle route problem to refine the earth allocation problem based on the genetic algorithm. In reality, earth allocation may include multiple vehicle replenishment centers (depots), multiple cutting fields (distribution centers), and multiple filling fields (customers) distributed geographically in a certain range of areas. Before reaching its mileage restriction, each vehicle must return to any replenishment site to refuel. Furthermore, fuel usage varies depending on whether the truck is empty or loaded. Furthermore, each excavation area has a maximum amount of earthwork it can generate. Hence, an excavation area cannot contribute a limitless amount of earthwork to its neighboring filling regions. When the earthwork from their nearest excavation location is depleted, those filling regions must choose another suitable excavation place to furnish the earthwork. In light of the aforementioned restrictions, the goal of this work is to maximize the number of sent trucks and the overall fuel usage throughout the shipping distance. Accordingly, this paper made the following contributions:(1)
This paper constructs a complicated multidepot two-echelon vehicle routing problem (MD-TEVRP) and provides a many-to-many recursive pairing solution based on the genetic algorithm.(2)
A self-adaptive mechanism is designed to control the crossover and mutation rate to manipulate and maintain the diversity of the generated population, which can prevent the local convergence problem and provide robust performance.(3)
The proposed method demonstrates great performance in the case study and comparative study, which provided guidance for the construction of the skiing courses of the Beijing Winter Olympic Games Skiing Center in Yanqing, Beijing.


## 2. Model Formation

### 2.1. Problem Description

This paper assumes that the earthwork working sites are divided into 4 parts, i.e., the replenishment centers, the cutting fields, the filling fields, and the paths. The replenishment centers, denoted asj∈J, are defined as the depots where engineering vehicles (in this paper, dump trucks are mainly considered) can get refueled or repaired. Each replenishment center is able to accommodate Ij dump trucks that can provide stable and nonstop hauling services at each dispatch, which means these dump trucks will not stop transporting earthworks from the cutting fields to the filling fields until they exhaust all fuel capacity. To this point, earthwork allocation can be regarded as an extension of the vehicle route problem (VRP) with the purpose of transporting the earthwork from several specific places to other different appointed locations. The planning key of earthwork allocation is to transport the required amount of earth volume from the cutting field k∈K and distribute them according to the different earthwork needs of the filling fields g∈G. To simplify this problem, this paper defines the volume of earthwork to be excavated in the kth cutting field as the number of fully loaded dump trucks bk needed to transport away all the dredged earth. Similarly, the volume of earthwork needed in the gth filling fields demands cg fully loaded dump trucks to unload the earth. It is notable that the earthwork from all cutting fields should be not less than the earthwork needed in all filling fields, namely,(1)∑k∈Kbk≥∑g∈Gcg,∀k∈K,∀g∈G,where bk is the total number of fully loaded dump trucks needed to haul all earthwork in kth cutting field, and cg is the total number of fully loaded dump trucks needed to fill the gth filling field.As illustrated in Figure2, an empty-laden truck is dispatched from the jth replenishment centers to the kth cutting field after fully loading with earth and rocks. This truck will move to the gth filling field and unload the earth and rocks. Depending on its fuel capacity, this truck can either return to one of the replenishment centers to get refueled (in this case, this path for this truck is called one-way path) or head to one of the cutting fields to get loaded again and repeat the earthwork hauling work between the cutting fields and the filling fields (in this case, the path for this truck is called multiway path). In this paper, the maximum mileage of a truck is dynamic based on how often it is fully loaded during each dispatch, i.e., each truck has its own maximum fuel capacity lij, and thus, a fully loaded truck will consume more fuel at each unit hauling distance, and the more frequent it is fully loaded, the less mileage it can guarantee. Practically, earthwork allocation is repetitive, and thus, a multiway transportation path planning for each dump truck is indispensable. On the basis of not exceeding the mileage limit, it is required that a minimum number of dispatched trucks and a global minimum hauling distance are derived accordingly to achieve a minimum fuel consumption for the least carbon emission.Figure 2
Schematic diagram of earthwork allocation (the dash lines indicate that the trucks have two options depending on their mileage limitation. They can either return to a replenishment center or head to a cutting field to reload).


### 2.2. Objective Model and Constraints

According to the problem illustrated above and the schematic diagram shown in Figure2, our objective is to find the shortest overall hauling distance and achieve the lowest fuel consumption. Therefore, the objective model can be formulated as follows:(2)minZ=∑j∈J∑i∈Ij∑k∈K∑g∈GαdijkXijjk+digkUijgk+digjVijgj+βdikgYijkg,where dijk is the hauling distance of the ith dump truck from the jth replenishment center to the kth cutting field, dikg is the hauling distance from the kth cutting field to the gth filling field, digk is the hauling distance from the gth filling field to the kth cutting field, and digj is the hauling distance from the gth filling field to the jth replenishment center. α and β are the average fuel consumption per unit distance of each dump truck at no load and full load, respectively. Xijjk, Yijkg, Uijgk, and Vijgj are binary decision variables.Equation (2) is the objective function that consists of two parts. The first part is the total fuel consumption of all dump trucks at no load. To be specific, αdijkXijjk indicates the fuel consumption of the ith dump truck at the jth replenishment center dispatched from the jth replenishment center to the kth cutting field, αdigkUijgk is the fuel consumption upon traveling from the gth filling field back to the kth cutting field at no load, and αdigjVijgj represents the fuel consumption coming from the gth filling field back to the jth replenishment center at no load. The second part is the total fuel consumption of all dump trucks at full load, specifically referring to the sum fuel consumption hauling from the kth cutting field to the gth filling field.In addition, the objective function in (1) should be subjected to the following constraints:(3)∑i∈Ij∑k∈KXijjk≤aj,∀i∈Ij,∀j∈J,∀k∈K,(4)∑j∈J∑i∈Ij∑g∈GXijjk+Uijgk≥bk,∀i∈Ij,∀j∈J,∀k∈K,∀g∈G,(5)∑j∈J∑i∈Ij∑k∈KYijkg≥cg,∀i∈Ij,∀j∈J,∀k∈K,∀g∈G,(6)∑k∈K∑g∈Gdijk+digk+dikjα+dikgβ≤lij,∀i∈Ij,∀j∈J,∀k∈K,∀g∈G,(7)∑k∈Kbk≥∑g∈Gcg,∀k∈K,∀g∈G,(8)Xijjk∈0,1,∀i∈Ij,∀j∈J,∀k∈K,(9)Yijkg∈0,1,∀i∈Ij,∀j∈J,∀k∈K,∀g∈G,(10)Uijgk∈0,1,∀i∈Ij,∀j∈J,∀k∈K,∀g∈G,(11)Vijgj∈0,1,∀i∈Ij,∀j∈J,∀g∈G,where aj is the total number of available dump trucks at the jth replenishment center. (3) represents that the dispatched dump trucks from each replenishment center should not exceed the maximum number of available dump trucks parked in each replenishment center. Equations (4) and (5) determine that the number of dump trucks sent to each cutting field and each filling field should satisfy the earthwork hauling need. Equation (6) indicates that the maximum hauling distance of each dump truck is constrained by its own maximum fuel capacity. Equation (7) indicates that the earth volume in all cutting fields must not be less than that in all filling fields, otherwise, the filling need will not be satisfied. Equations (8) to (11) are the decision variables to control the hauling path of each dump truck.


## 2.1. Problem Description

This paper assumes that the earthwork working sites are divided into 4 parts, i.e., the replenishment centers, the cutting fields, the filling fields, and the paths. The replenishment centers, denoted asj∈J, are defined as the depots where engineering vehicles (in this paper, dump trucks are mainly considered) can get refueled or repaired. Each replenishment center is able to accommodate Ij dump trucks that can provide stable and nonstop hauling services at each dispatch, which means these dump trucks will not stop transporting earthworks from the cutting fields to the filling fields until they exhaust all fuel capacity. To this point, earthwork allocation can be regarded as an extension of the vehicle route problem (VRP) with the purpose of transporting the earthwork from several specific places to other different appointed locations. The planning key of earthwork allocation is to transport the required amount of earth volume from the cutting field k∈K and distribute them according to the different earthwork needs of the filling fields g∈G. To simplify this problem, this paper defines the volume of earthwork to be excavated in the kth cutting field as the number of fully loaded dump trucks bk needed to transport away all the dredged earth. Similarly, the volume of earthwork needed in the gth filling fields demands cg fully loaded dump trucks to unload the earth. It is notable that the earthwork from all cutting fields should be not less than the earthwork needed in all filling fields, namely,(1)∑k∈Kbk≥∑g∈Gcg,∀k∈K,∀g∈G,where bk is the total number of fully loaded dump trucks needed to haul all earthwork in kth cutting field, and cg is the total number of fully loaded dump trucks needed to fill the gth filling field.As illustrated in Figure2, an empty-laden truck is dispatched from the jth replenishment centers to the kth cutting field after fully loading with earth and rocks. This truck will move to the gth filling field and unload the earth and rocks. Depending on its fuel capacity, this truck can either return to one of the replenishment centers to get refueled (in this case, this path for this truck is called one-way path) or head to one of the cutting fields to get loaded again and repeat the earthwork hauling work between the cutting fields and the filling fields (in this case, the path for this truck is called multiway path). In this paper, the maximum mileage of a truck is dynamic based on how often it is fully loaded during each dispatch, i.e., each truck has its own maximum fuel capacity lij, and thus, a fully loaded truck will consume more fuel at each unit hauling distance, and the more frequent it is fully loaded, the less mileage it can guarantee. Practically, earthwork allocation is repetitive, and thus, a multiway transportation path planning for each dump truck is indispensable. On the basis of not exceeding the mileage limit, it is required that a minimum number of dispatched trucks and a global minimum hauling distance are derived accordingly to achieve a minimum fuel consumption for the least carbon emission.Figure 2
Schematic diagram of earthwork allocation (the dash lines indicate that the trucks have two options depending on their mileage limitation. They can either return to a replenishment center or head to a cutting field to reload).


## 2.2. Objective Model and Constraints

According to the problem illustrated above and the schematic diagram shown in Figure2, our objective is to find the shortest overall hauling distance and achieve the lowest fuel consumption. Therefore, the objective model can be formulated as follows:(2)minZ=∑j∈J∑i∈Ij∑k∈K∑g∈GαdijkXijjk+digkUijgk+digjVijgj+βdikgYijkg,where dijk is the hauling distance of the ith dump truck from the jth replenishment center to the kth cutting field, dikg is the hauling distance from the kth cutting field to the gth filling field, digk is the hauling distance from the gth filling field to the kth cutting field, and digj is the hauling distance from the gth filling field to the jth replenishment center. α and β are the average fuel consumption per unit distance of each dump truck at no load and full load, respectively. Xijjk, Yijkg, Uijgk, and Vijgj are binary decision variables.Equation (2) is the objective function that consists of two parts. The first part is the total fuel consumption of all dump trucks at no load. To be specific, αdijkXijjk indicates the fuel consumption of the ith dump truck at the jth replenishment center dispatched from the jth replenishment center to the kth cutting field, αdigkUijgk is the fuel consumption upon traveling from the gth filling field back to the kth cutting field at no load, and αdigjVijgj represents the fuel consumption coming from the gth filling field back to the jth replenishment center at no load. The second part is the total fuel consumption of all dump trucks at full load, specifically referring to the sum fuel consumption hauling from the kth cutting field to the gth filling field.In addition, the objective function in (1) should be subjected to the following constraints:(3)∑i∈Ij∑k∈KXijjk≤aj,∀i∈Ij,∀j∈J,∀k∈K,(4)∑j∈J∑i∈Ij∑g∈GXijjk+Uijgk≥bk,∀i∈Ij,∀j∈J,∀k∈K,∀g∈G,(5)∑j∈J∑i∈Ij∑k∈KYijkg≥cg,∀i∈Ij,∀j∈J,∀k∈K,∀g∈G,(6)∑k∈K∑g∈Gdijk+digk+dikjα+dikgβ≤lij,∀i∈Ij,∀j∈J,∀k∈K,∀g∈G,(7)∑k∈Kbk≥∑g∈Gcg,∀k∈K,∀g∈G,(8)Xijjk∈0,1,∀i∈Ij,∀j∈J,∀k∈K,(9)Yijkg∈0,1,∀i∈Ij,∀j∈J,∀k∈K,∀g∈G,(10)Uijgk∈0,1,∀i∈Ij,∀j∈J,∀k∈K,∀g∈G,(11)Vijgj∈0,1,∀i∈Ij,∀j∈J,∀g∈G,where aj is the total number of available dump trucks at the jth replenishment center. (3) represents that the dispatched dump trucks from each replenishment center should not exceed the maximum number of available dump trucks parked in each replenishment center. Equations (4) and (5) determine that the number of dump trucks sent to each cutting field and each filling field should satisfy the earthwork hauling need. Equation (6) indicates that the maximum hauling distance of each dump truck is constrained by its own maximum fuel capacity. Equation (7) indicates that the earth volume in all cutting fields must not be less than that in all filling fields, otherwise, the filling need will not be satisfied. Equations (8) to (11) are the decision variables to control the hauling path of each dump truck.


## 3. Adaptive Genetic Algorithm

The proposed earthwork allocation problem can be regarded as an MD-TEVRP, which is a many (dump trucks in the replenishment centers)-to-many (filling fields) matching problem, considering the transfer stations (cutting fields). The MD-TEVRP is an NP-hard question [41], and thus, this paper proposed an adaptive genetic algorithm to determine the optimal solution.The genetic algorithm is a heuristic algorithm that is based on the “Survival of the Fittest” and “Natural Selection” theories. Figure3 depicts the process and technique by which randomly created individuals compete with one another and form new generations using self-adaptive crossover and mutation strategies. When the system converges to a stable solution in the final phase, it signifies that the optimal solution has been discovered.Figure 3
Flow chart of the adaptive genetic algorithm.


### 3.1. Generation of Initial Population

Based on the dynamic planning theory, the Floyd algorithm is used to calculate the shortest distance between the arbitrary pointsi and j in the original spatial weight matrix [42]. The Floyd algorithm seeks a third point k between the points i and j by comparing the distances from i to j directly and from i to j via k. If the latter distance is shorter, the path i-k-j is updated as the shortest distance between the points i and j. Therefore, the recursive process of the Floyd algorithm can be expressed as follows:(12)dijk=mindijk−1,dikk−1+djkk−1.By deriving the shortest distance between each point via (4), the next step is to find the optimal path for the many-to-many problem described in Section 2. As depicted in Figure 3, the next step is to determine the coding of the initial population.Assume that there areN units (chromosomes) in each generation, and each chromosome contains two parts: one is the feasible solution, while the other is the fitness value. Therefore, the information of a chromosome can be stored in a 1×2 Cell Array. According to the proposed model, the coding for the first element in the array can be designed as a ∑g∈Gcg×7 matrix to store the information of a feasible solution. The formation of this matrix is generated by our proposed many-to-many algorithm, as seen in Table 1.Table 1
Algorithm for generating initial population.
#pseudocodeLoad node information;Load distance informationD and path information P from Floyd algorithm;for each unitn∈Nfor each cutting fieldk∈Kifbk−∑j∈J∑i∈Ij∑g∈GXijjk+Uijgk≥0ifaj−∑i∈Ij∑k∈KXijjk≥0Randomly select a replenishment centerj∈J;end ififcg−∑j∈J∑i∈Ij∑k∈KYjkg≥0Randomly select a filling fieldg∈G;End ifEnd forConstruct a one-way path solutionS as shown in Table 1;While∑k∈K∑g∈Gdijk+digk+dikjα+dikgβ<lijS: ,6←0, S:,2←x and append1row⟶S;End whileEnd forTo be more specific, if the dump truck was assigned a one-way path because of its mileage limit, the one-way path solution can be expressed as a1×7 matrix that contains the following information (see Table 2):Table 2
An example of a one-way path feasible solution.
1024681Indexing the number of dispatched dump trucksIndexing it is a one-way path solutionThe serial number of the departure replenishment centerThe serial number of the cutting field arrived atThe serial number of the filling field reachedThe serial number of the returned replenishment centerThe total distance the truck travels in the current dispatchSimilarly, a multiway path viable solution can be expressed in ax×7 matrix shown in Table 3. It is noticeable that the final retuning replenishment center should be the same as the departure one.Table 3
An example of a multiway path feasible solution.
Indexing the number of dispatched dump trucksIndexing how many runs the current truck travelsThe serial number of the departure replenishment centerThe serial number of the cutting field arrived atThe serial number of the filling field reached0 if the dump trucks can continue shipment, otherwise indexing the serial number of the returned replenishment centerThe distance the truck travels in each run112560512479010…1x325214212560522479010…2x325214…The second element of the array records the fitness value of each chromosome, which stores the sum of the total fuel consumption of all dispatched trucks in the corresponding feasible solution.


### 3.2. Roulette Select

Roulette selection is to determine some better individuals from the paternal chromosome based on the fitness value associated with the paternal chromosome. The smaller the fitness value, the more likely the corresponding chromosome is to be inherited. In addition, the elite strategy was adopted for the offspring generated after the above random operation, in which the best feasible solution of the parent gene was retained and the worst feasible solution of the offspring was replaced by the best feasible solution of the parent. The fitness value is defined as follows:(13)fn=1∑j∈J∑i∈Ij∑k∈K∑g∈GαdijkXijjk+digkUijgk+digjVijgj+βdikgYijkg,n∈N.The probability of selecting each unit from the population is calculated as(14)pn=fn∑n∈Nfn.


### 3.3. Crossover

In this model, the individual pairs of chromosomes can intersect with each other according to the crossover rate defined in (8), i.e., a certain rate of individual pairs can exchange their cutting fields and recalculate their overall mileage and fitness value thereafter. The crossover rate is self-adaptive according to the proportion p of the ratio of current optimal individuals to the whole population size. ω is the allowable peak crossover rate, and thus, CR∈0,ω∀p∈0,1.(15)CR=ωlog5050p.Specifically, if the crossover rate is satisfied, thenth unit will swap its cutting field with the (n + N/2)th unit. The algorithm is presented in Table 4.Table 4
Algorithm for crossover.
#pseudocodeFor each unitn∈N/2ifrandom0,1≤CRSn,4⇄Sn+N/2,4;end ifRecomputefn=1/∑j∈J∑i∈Ij∑k∈K∑g∈GαdijkXijjk+digkUijgk+digjVijgj+βdikgYijkg;end forEquation (8) is an increasing logarithmic function subjected to the proportion of the current optimal individuals. The shape of (8) is demonstrated in Figure 4. As indicated, as the number of optimal individual increases, the crossover rate will rise to increase the population diversity to avoid premature problem. It is noticeable that when p is larger than 0.02, the crossover rate and the mutation rate will become positive. p is the proportion of the current optimal individuals to the whole population. In this case, we assume that only when the proportion of optimal individuals is larger than 2%, the self-adaptive mechanism will work and start to manipulate the diversity of population. Otherwise, the algorithm will let the nature take its course.Figure 4
Logarithmic functions with a base of 50.


### 3.4. Mutation

In the mutation stage, this paper mainly carries out the variation operations on the fourth element (cutting fields) of the first array. Specifically, this paper will replace the cutting field with another available one with the remaining earth volume required to be transported, which can achieve a shorter hauling distance from the replenishment centers to the filling fields on this single run. Therefore, it is necessary to ensure that the amount of earthwork in the cutting fields is not less than that in the filling fields when modeling, so that it can provide redundant earth volume in the cutting fields for variation. After mutation, the path length traveled by the dump trucks and the corresponding fitness value of the entire feasible solution will be updated again. The mutation rate is also controlled by a self-adaptive mechanism as (16), where χ is the maximum mutation rate, which ensures MR∈0,χ∀p∈0,1.(16)MR=χlog5050p.As displayed in Table5, the mutation algorithm aims to find the best cutting field to achieve a lower fitness value for each run. Once the mutation rate is satisfied, the algorithm will find all available cutting fields that can replace the original one and select the one with the least fitness value among all the alternatives.Table 5
Algorithm for mutation.
#pseudocodeFor each unitn∈Nifrandom0,1≤MRKremain=find available cutting fields with remaining earth volume inK;forkremain∈KremainSn,4←kremain;Recomputefn=1/∑j∈J∑i∈Ij∑k∈K∑g∈GαdijkXijjk+digkUijgk+digjVijgj+βdikgYijkg;end forfnmink=minfnkremain;Sn,4←mink;end ifend forThe adaptive mechanism in the mutation stage is similar to that in the crossover stage, which aims to enlarge the diversity of the population and prevent a premature problem. As shown in Figure5, during 4500 iterations, the optimal solution at each generation oscillates up and down greatly, ranging from 31,500 to 33,750, which means a high divergence population is maintained by this self-adaptive mechanism to avoid a premature problem. This adaptive mechanism provides more possibility to seek better solutions by jumping out local convergence.Figure 5
High diversity of optimal solutions at each generation.


### 3.5. Recursive Optimization

A better solution with the lower fitness value can be obtained by repeating sections3.2 to 3.4 and iteratively updating the chromosomal information of each generation of the population until meeting the following convergence criterion:(17)number of iterations=max500,ns,where ns is the number of iterations when the optimal solution remains unchanged for over 50 iterations.


## 3.1. Generation of Initial Population

Based on the dynamic planning theory, the Floyd algorithm is used to calculate the shortest distance between the arbitrary pointsi and j in the original spatial weight matrix [42]. The Floyd algorithm seeks a third point k between the points i and j by comparing the distances from i to j directly and from i to j via k. If the latter distance is shorter, the path i-k-j is updated as the shortest distance between the points i and j. Therefore, the recursive process of the Floyd algorithm can be expressed as follows:(12)dijk=mindijk−1,dikk−1+djkk−1.By deriving the shortest distance between each point via (4), the next step is to find the optimal path for the many-to-many problem described in Section 2. As depicted in Figure 3, the next step is to determine the coding of the initial population.Assume that there areN units (chromosomes) in each generation, and each chromosome contains two parts: one is the feasible solution, while the other is the fitness value. Therefore, the information of a chromosome can be stored in a 1×2 Cell Array. According to the proposed model, the coding for the first element in the array can be designed as a ∑g∈Gcg×7 matrix to store the information of a feasible solution. The formation of this matrix is generated by our proposed many-to-many algorithm, as seen in Table 1.Table 1
Algorithm for generating initial population.
#pseudocodeLoad node information;Load distance informationD and path information P from Floyd algorithm;for each unitn∈Nfor each cutting fieldk∈Kifbk−∑j∈J∑i∈Ij∑g∈GXijjk+Uijgk≥0ifaj−∑i∈Ij∑k∈KXijjk≥0Randomly select a replenishment centerj∈J;end ififcg−∑j∈J∑i∈Ij∑k∈KYjkg≥0Randomly select a filling fieldg∈G;End ifEnd forConstruct a one-way path solutionS as shown in Table 1;While∑k∈K∑g∈Gdijk+digk+dikjα+dikgβ<lijS: ,6←0, S:,2←x and append1row⟶S;End whileEnd forTo be more specific, if the dump truck was assigned a one-way path because of its mileage limit, the one-way path solution can be expressed as a1×7 matrix that contains the following information (see Table 2):Table 2
An example of a one-way path feasible solution.
1024681Indexing the number of dispatched dump trucksIndexing it is a one-way path solutionThe serial number of the departure replenishment centerThe serial number of the cutting field arrived atThe serial number of the filling field reachedThe serial number of the returned replenishment centerThe total distance the truck travels in the current dispatchSimilarly, a multiway path viable solution can be expressed in ax×7 matrix shown in Table 3. It is noticeable that the final retuning replenishment center should be the same as the departure one.Table 3
An example of a multiway path feasible solution.
Indexing the number of dispatched dump trucksIndexing how many runs the current truck travelsThe serial number of the departure replenishment centerThe serial number of the cutting field arrived atThe serial number of the filling field reached0 if the dump trucks can continue shipment, otherwise indexing the serial number of the returned replenishment centerThe distance the truck travels in each run112560512479010…1x325214212560522479010…2x325214…The second element of the array records the fitness value of each chromosome, which stores the sum of the total fuel consumption of all dispatched trucks in the corresponding feasible solution.


## 3.2. Roulette Select

Roulette selection is to determine some better individuals from the paternal chromosome based on the fitness value associated with the paternal chromosome. The smaller the fitness value, the more likely the corresponding chromosome is to be inherited. In addition, the elite strategy was adopted for the offspring generated after the above random operation, in which the best feasible solution of the parent gene was retained and the worst feasible solution of the offspring was replaced by the best feasible solution of the parent. The fitness value is defined as follows:(13)fn=1∑j∈J∑i∈Ij∑k∈K∑g∈GαdijkXijjk+digkUijgk+digjVijgj+βdikgYijkg,n∈N.The probability of selecting each unit from the population is calculated as(14)pn=fn∑n∈Nfn.


## 3.3. Crossover

In this model, the individual pairs of chromosomes can intersect with each other according to the crossover rate defined in (8), i.e., a certain rate of individual pairs can exchange their cutting fields and recalculate their overall mileage and fitness value thereafter. The crossover rate is self-adaptive according to the proportion p of the ratio of current optimal individuals to the whole population size. ω is the allowable peak crossover rate, and thus, CR∈0,ω∀p∈0,1.(15)CR=ωlog5050p.Specifically, if the crossover rate is satisfied, thenth unit will swap its cutting field with the (n + N/2)th unit. The algorithm is presented in Table 4.Table 4
Algorithm for crossover.
#pseudocodeFor each unitn∈N/2ifrandom0,1≤CRSn,4⇄Sn+N/2,4;end ifRecomputefn=1/∑j∈J∑i∈Ij∑k∈K∑g∈GαdijkXijjk+digkUijgk+digjVijgj+βdikgYijkg;end forEquation (8) is an increasing logarithmic function subjected to the proportion of the current optimal individuals. The shape of (8) is demonstrated in Figure 4. As indicated, as the number of optimal individual increases, the crossover rate will rise to increase the population diversity to avoid premature problem. It is noticeable that when p is larger than 0.02, the crossover rate and the mutation rate will become positive. p is the proportion of the current optimal individuals to the whole population. In this case, we assume that only when the proportion of optimal individuals is larger than 2%, the self-adaptive mechanism will work and start to manipulate the diversity of population. Otherwise, the algorithm will let the nature take its course.Figure 4
Logarithmic functions with a base of 50.


## 3.4. Mutation

In the mutation stage, this paper mainly carries out the variation operations on the fourth element (cutting fields) of the first array. Specifically, this paper will replace the cutting field with another available one with the remaining earth volume required to be transported, which can achieve a shorter hauling distance from the replenishment centers to the filling fields on this single run. Therefore, it is necessary to ensure that the amount of earthwork in the cutting fields is not less than that in the filling fields when modeling, so that it can provide redundant earth volume in the cutting fields for variation. After mutation, the path length traveled by the dump trucks and the corresponding fitness value of the entire feasible solution will be updated again. The mutation rate is also controlled by a self-adaptive mechanism as (16), where χ is the maximum mutation rate, which ensures MR∈0,χ∀p∈0,1.(16)MR=χlog5050p.As displayed in Table5, the mutation algorithm aims to find the best cutting field to achieve a lower fitness value for each run. Once the mutation rate is satisfied, the algorithm will find all available cutting fields that can replace the original one and select the one with the least fitness value among all the alternatives.Table 5
Algorithm for mutation.
#pseudocodeFor each unitn∈Nifrandom0,1≤MRKremain=find available cutting fields with remaining earth volume inK;forkremain∈KremainSn,4←kremain;Recomputefn=1/∑j∈J∑i∈Ij∑k∈K∑g∈GαdijkXijjk+digkUijgk+digjVijgj+βdikgYijkg;end forfnmink=minfnkremain;Sn,4←mink;end ifend forThe adaptive mechanism in the mutation stage is similar to that in the crossover stage, which aims to enlarge the diversity of the population and prevent a premature problem. As shown in Figure5, during 4500 iterations, the optimal solution at each generation oscillates up and down greatly, ranging from 31,500 to 33,750, which means a high divergence population is maintained by this self-adaptive mechanism to avoid a premature problem. This adaptive mechanism provides more possibility to seek better solutions by jumping out local convergence.Figure 5
High diversity of optimal solutions at each generation.


## 3.5. Recursive Optimization

A better solution with the lower fitness value can be obtained by repeating sections3.2 to 3.4 and iteratively updating the chromosomal information of each generation of the population until meeting the following convergence criterion:(17)number of iterations=max500,ns,where ns is the number of iterations when the optimal solution remains unchanged for over 50 iterations.


## 4. Validation of the Proposed Algorithm

To validate the proposed algorithm, an earthwork allocation project on the G1 and D2 skiing courses of the Beijing Winter Olympic Games Skiing Center in 2022 is used. The schematic diagram of the original transportation network is shown in Figure6, and the initial conditions of this case study are listed as follows:(1)
Overall excavation volume: 2352 m3(2)
Overall filling volume: 1978 m3(3)
Maximum carrying capacity of each dump truck: 10 m3(4)
The transportation fuel consumption of empty load is 0.1 L/km(5)
The transportation fuel consumption of full loaded is 0.2 L/km(6)
The maximum capacity of fuel tank: 100 L(7)
Number of cutting fields: 3(8)
Number of filling fields: 7(9)
Number of replenishment centers: 2Figure 6
Skiing courses of the Beijing Winter Olympic Games Skiing Center in Yanqing District. (a) Real map of the skiing courses. (b) Distribution plan of the skiing courses. (c) BIM map of the G1 and D2 skiing courses. (d) Schematic diagram of the transportation network between the two depots, G1 and D2 skiing courses.
(a)(b)(c)(d)The node information is listed in Table6. The number of dump trucks in Table 6 provides the meanings, subject to the locations, respectively, as follows:(1)
Number of dump trucks each replenishment center can provide(2)
Number of dump trucks needed to transport all earth volume away from the cutting fields(3)
Number of dump trucks needed to transport the required earth volume to the filling fieldsTable 6
Node information.
Replenishment centersCutting fieldsFilling fieldsNode123456789101112Number of vehicles4030835785732280260330300250334224According to the distance information read from Figure3, the spatial weight matrix can be constructed as Table 7.Table 7
The spatial weight matrix.
DestinationsOrigins0—36—2368—11—0—0—231———3—5—3—0—---—1——4—362—0—62———946—3——0————6—2—21—6—0———17—23——2——0—6——336—1————02———68—————620———8—3——61———0———11—49—7————0—11—5—42—3————0—0—36—2368—11—0According to equation (12), the distance information matrix and the path information matrix can be derived (see Tables 8 and 9).Table 8
The Floyd distance information matrix (1 = 5 km).
NodeDestinationsNode123456789101112Origins103356234637623062314792853360895613649452805329849456395045101151126215340568176734625507661038471910670275109693811862097910326451679087117849117105780131265942631097130Table 9
The Floyd path information matrix.
NodeDestinationsNode123456789101112Origins11636667336372626456466661231131111881111422242272721112522225212212221261212261111011271414121719111283333333893339888778789887106666666661066113634663336114127274527772412The overall demand of all filling areas is 1978, and thus, the final optimal solution is a1978×7 matrix. Each row represents the path information of a single run. Table 10 demonstrates the results of the optimal solution for the last 20 rows. It can be seen that 2 dump trucks are dispatched for multiple deliveries (e.g., the 63rd dump truck has been assigned to carry out 12 transportation tasks from 1959 to 1970). Therefore, the optimal path can be read assisted by Table 9 generated from Floyd Algorithm. For instance, the 1959th run travels from 1-4-11 as indicated in Table 10; however, the real hauling path from 1 to 11 is 1-6-2-4-11 by reading the path information matrix in Table 9. It is worth noting that 64 dump trucks are used to carry out the earthwork allocation task and achieve minimal fuel consumption so far.Table 10
Optimal solution of the last 20 rows.
Indexing the number of dispatched dump trucksIndexing how many runs of current truck travelsThe serial number of the departure replenishment centerThe serial number of the cutting field arrived atThe serial number of the filling field reached0 if the dump trucks can continue shipment, otherwise indexing the serial number of the returned replenishment centerThe distance the truck travels in each run63195914110146319601149017631961949016631962949016631963949016631964949016631965949016631966949016631967949016631968949016631969949016631970949122641971149013641972959022641973959022641974959022641975959022641976959022641977959022641978949122Meanwhile, Figure7 is the iterative optimizing process. The iteration operates 500 times in total and converges at the 440th generation. The final optimal fuel consumption is 31352 L.Figure 7
Recursive optimizing process of earthwork allocation.To further validate the performance of our proposed algorithm, a comparative study has been made by comparing our method with former work in [6]. The result in Figure 8 demonstrated that our method showed a dominant advantage over the past work, and roughly, a performance boost of 38% is achieved. The global optimal hauling cost in Figure 8 is related to fuel consumption. Thus, a slight change to the algorithm is made to accord with the optimization objective.Figure 8
Comparative study between [6] and our method. (a) The former work with two-way path planning. (b) Our method with multiway path planning.
(a)(b)


## 5. Conclusions

This paper constructs a complex earthwork allocation path-planning modeling based on MD-TEVRP by considering the multiple runs of dump trucks between the cutting fields and the filling fields. A well-designed genetic algorithm based on this model is proposed to achieve the following contributions:(1)
This paper provides a multiway transportation solution to a typical MD-TEVRP, which achieves many-to-many recursive pairing in earthwork allocation path planning.(2)
The proposed self-adaptive mechanism can manipulate the diversity of the generated populations to prevent the premature problem of GA and achieve a lower fitness value.(3)
The optimization method shows a dominant advantage over the past work by increasing approximately 38% of performance. The proposed method also provides guidance for the construction design of the Skiing Courses of Beijing Winter Olympic Games Skiing Center in Yanqing, Beijing.The implication of variables is listed in Table11.Table 11
Variables and explanation.
VariablesExplanationJNumber of replenishment centersIjNumber of dump trucks in thejth replenishment centerKNumber of cutting fieldsGNumber of filling fieldsdijkThe hauling distance of theith dump truck from the jth replenishment center to the kth cutting fielddikgThe hauling distance of theith dump truck from the kth cutting field to the gth filling fielddigkThe hauling distance of theith dump truck from the gth filling field to the kth cutting fielddigjThe hauling distance of theith dump truck from the gth filling field to the jth replenishment centerαThe average fuel consumption per unit distance of each dump truck at no loadβThe average fuel consumption per unit distance of each dump truck at full loadXijjkBinary decision variable: whether theith dump truck at the jth replenishment center travels from the jth replenishment center to the kth cutting fieldYijkgBinary decision variable: whether theith dump truck at the jth replenishment center travels from the kth cutting field to the gth filling fieldUijgkBinary decision variable: whether theith dump truck at the jth replenishment center travels from the gth filling field to the kth cutting fieldVijgjBinary decision variable: whether theith dump truck at the jth replenishment center travels from the gth filling field to the jth replenishment centerajTotal number of available dump trucks at thejth replenishment centerbkTotal number of fully loaded dump trucks needed to haul all earthwork in thekth cutting fieldcgTotal number of fully loaded dump trucks needed to fill thegth filling fieldlijThe maximum fuel capacity of theith dump truck at the jth replenishment centerNPopulation sizefnThe fitness value of the nth unitpnThe probability of selecting each unit from the populationωThe allowable peak crossover rateCRThe crossover ratepThe ratio of current optimal individuals to the whole population sizeχThe maximum mutation rateMRThe mutation ratensThe number of iterations when the optimal solution remains unchanged for over 50 iterations



---

*Source: 8373138-2022-01-29.xml*