Time window constrained vehicle routing method based on multi-objective evolutionary algorithm

Abstract

The application discloses a time window constraint vehicle path planning method based on a multi-target evolution algorithm, which comprises a balance separation control mechanism and a Kanpo two-stage strategy model bearing the mechanism; an initial population is obtained by adopting a random traversal insertion initialization method, and a search tendency is adaptively adjusted by the established Kanpo two-stage strategy model; a mating pool is generated by a stage-type mating selection operation matching the characteristics of the search tendency, an ordered crossover and a reversal mutation operator are adopted, and a candidate offspring set is generated; a next generation population is screened out by a stage-type environment selection operation matching the characteristics of the search tendency; and an optimal solution set generated by evolution can obtain a transport scheme with extremely high utility performance, and can effectively solve the problem of invalid decision scheme caused by conflicts between targets under a complex background.

Description

TECHNICAL FIELD

[0001] The application belongs to the technical field of multi-objective optimization research in evolutionary intelligent evolutionary computation, and particularly relates to a time window constrained vehicle path planning method based on a multi-objective evolutionary algorithm. BACKGROUND

[0002] Vehicle path planning has a complex problem background, and specifically relates to multiple constraint conditions and multiple target functions that need to be optimized. On the premise of not exceeding the vehicle loading capacity, the decision maker expects the transportation scheme to meet the customer's loading demand to the greatest extent with the minimum comprehensive cost. However, the conflict between the targets makes optimization difficult, and fewer vehicle numbers often mean smaller total loading capacity.

[0003] The multi-objective evolutionary algorithm is a heuristic numerical random search algorithm based on Darwin's theory of evolution. This kind of algorithm can use problem information to guide the search to obtain the numerical optimal solution. However, as one of the important methods for solving the vehicle path planning problem, the multi-objective evolutionary algorithm presents a sharp decline in utility performance when dealing with problems involving three or more target functions. The decision maker may use a transportation scheme with low utility performance, which may lead to problems such as accelerated aging of some vehicles, higher comprehensive cost than actually needed, serious violation of customer loading and unloading requirements and service time window, and thus adversely affect the service provider in decision management, equipment maintenance, fund balance, employee and customer satisfaction, etc. The reason for the above problems is that the compromise search tendency of the existing algorithm is difficult to balance the conflict between the targets, and the mutual influence between different search tendencies is difficult to be completely separated, thus seriously deteriorating the decision result of the transportation scheme in a certain target. SUMMARY

[0004] The purpose of the present application is to provide a time window constrained vehicle path planning method based on a multi-objective evolutionary algorithm, which solves the problem of invalid decision scheme in the complex problem background in the prior art.

[0005] The technical solution adopted by the present application is: a time window constrained vehicle path planning method based on a multi-objective evolutionary algorithm, which is implemented according to the following steps:

[0006] Step 1: an initial population is obtained by using a random traversal insertion initialization method, the initial population represents an initial solution set of the time window constrained vehicle path planning problem, and the chromosome of any solution in the initial solution set represents the decision vector of the solution;

[0007] Step 2: a reference vector set is generated by using a hyperplane construction method, the reference vector set contains the same number of reference vectors as the number of solutions in the initial population, and is used to assist in the calculation of state representation values and the implementation of control mechanisms;

[0008] Step 3, the initial population is subjected to cyclic iterative evolution, and the evolution stage of the current population is obtained by the established two-stage strategy model of Campanelli, and the search tendency is adaptively adjusted, and the timing of determining the search tendency switching is configured to guarantee the effectiveness of the balanced separation;

[0009] Step 4, the evolution stage obtained in step 3 is used to select the dominance relationship and the characteristic value that match the characteristics of the stage to indicate the state of the solution set, and the mating pool is generated by the binary tournament with the convergence measure as the primary selection standard and the diversity measure as the secondary selection standard;

[0010] Step 5, the candidate offspring set is generated by using the ordered crossover and inversion mutation operators on the mating pool obtained in step 4;

[0011] Step 6, the evolution stage obtained in step 3 is used to select the dominance relationship and the characteristic value that match the characteristics of the stage to indicate the state of the solution set, and the next generation population is selected from the union set of the current population and the candidate offspring set obtained in step 5 by using the matching maintenance strategy;

[0012] Step 7, it is judged whether the evolution number of the next generation population obtained in step 6 reaches the set maximum evolution number, if not, it is indicated that the termination condition is not met, and the process returns to step 3, otherwise, the next generation population obtained in step 6 is output as the optimal solution set.

[0013] The application is also characterized in that,

[0014] As a preferred technical solution of the application, the two-stage strategy model of Campanelli in step 3 includes a set of balanced separation control mechanisms, including first, second and third levels; in the first level, the balanced state of the solution set is separated into convergence state and distribution state, and the corresponding compromise search tendency is separated into convergence emphasis and diversity emphasis, and the stages with different search tendencies are alternately executed according to the switching rules of the automatic configuration component; in the second level, the stages with different search tendencies correspondingly use different emphasis of the dominance rules, selection standards and maintenance strategies in the dominance relationship, mating selection component and environment selection component; in the third level, the convergence state and distribution state of different stages are characterized by the measurement values that match the characteristics of the stage.

[0015] As a preferred technical solution of the application, the automatic configuration component of the two-stage strategy model of Campanelli in step 3 is based on statistical experiments; the switching rule is obtained by simulating the macroeconomic regulation and control method of the Keynesian school and analyzing the ratio parameter μ of the first stage evolution number and the maximum evolution number, and in at least the first 75% of the search process, the automatic configuration component guides the search to stay in the first stage which emphasizes convergence, and in the last 25% of the search process, the automatic configuration component guides the two-stage alternately according to the number of objective functions and certain randomness.

[0016] In step 3, the number of objective functions and the current generation are obtained, and the current evolutionary stage is calculated using formula (1):

[0017]

[0018] Where I represents the current evolutionary stage, G done G represents the current generation number. max M represents the maximum number of generations, and M represents the number of objective functions. max With M min represents the maximum and minimum values ​​of the objective function in the statistical experiment, respectively, and rand(1) represents a random number in the interval [0, 1].

[0019] As a preferred embodiment of the present invention, the dominance relationships in step 4 that match the stage characteristics are the hybrid sub-region dominance relationship and the dual-grid dominance relationship:

[0020] Regarding the dominance relation of mixed subregions: Assuming u and v are two solutions, u dominates v when one of the following conditions is true, denoted as .

[0021] 1)

[0022] 2) u and v are equivalent in the Pareto sense;

[0023] a) Under normalization, u and v are associated with the same reference vector;

[0024] 1.

[0025] 2. However, under standardization, u and v are related to the same reference vector; PBI(u) < PBI(v), d1(u) < d1(v) and

[0026] b) Under normalization, u and v are associated with different reference vectors, but under standardization, u and v are associated with the same reference vector; And PBI(u)<PBI(v), And d1(u) < d1(v),

[0027] in The label indicates the normalized metric value. The symbol indicates Pareto dominance; d1 and d2 represent the projected distance and perpendicular distance of the solution on its associated reference vector; SD represents the density of the subregion to which the solution belongs, i.e., the number of solutions associated with the same reference vector; PBI represents the metric obtained using the penalty-based boundary crossing method and formula (2) in the MOEA / D algorithm designed by Zhang et al. in 2007.

[0028] PBI(x)=d1(x)+θd2(x), θ=1 (2)

[0029] Regarding the dominance relation of a double grid: Assuming u and v are two solutions, u dominates v when one of the following conditions is true, denoted as .

[0030] 1)

[0031] 2) u and v are equivalent in the Pareto sense;

[0032] a) Under the first-level grid system

[0033] b) u and v are equivalent in the first-level grid system; GS(u) < GS(v) and Level(u) < Level(v);

[0034] in The label indicates the grid dominance used in the GrEA algorithm, i.e., the first-level grid system; GS represents the sum of the coordinate parameters of the hyperbox to which the solution belongs, and solutions within the same hyperbox share the GS metric value; Level represents the horizontal level of the hyperbox, and hyperboxes at the same horizontal level share the Level metric value.

[0035] As a preferred technical solution of the present invention, the dual mesh dominance relationship in step 4 achieves adaptive adjustment of the superbox size through formula (3):

[0036]

[0037] Where (||·||) represents the number of internal elements to be calculated, div represents the number of cuts to uniformly divide each objective axis, α is a threshold factor related to the number of objective functions M, α=2M(M-2), P t Let denot be the population of generation t, and PF1 be the set of non-dominated solutions in the population of generation t that are on the first Pareto front.

[0038] As a preferred technical solution of the present invention, step 4 specifically includes:

[0039] Step 4.1: Match the applicable dominance relation and representation value according to the evolutionary stage obtained in Step 3. If the evolution is in the first stage, select the mixed sub-region dominance relation and d2 distance metric; otherwise, select the dual grid dominance relation and SD density metric. Obtain the non-dominance level of the current solution set through the non-dominance sorting method based on the corresponding dominance relation and store it in the convergence metric CM. Store the metric value of the current solution set in the diversity metric DM.

[0040] Step 4.2: Randomly select solution a and solution b from the current solution set and compare them according to the convergence metric CM and diversity metric DM obtained in Step 4.1. The solution with the smaller CM value wins the comparison; if the CM values ​​are the same, the solution with the smaller DM value wins the comparison; otherwise, randomly select one of the two to win. Add the winning solution to the mating pool and repeat the above comparison operation until the mating pool is saturated.

[0041] As a preferred technical solution of the present invention, step 6 specifically includes:

[0042] Step 6.1: Match the applicable dominance relation, representation value, and maintenance strategy according to the evolutionary stage obtained in Step 3. If the evolution is in the first stage, select the mixed sub-region dominance relation and d2 distance metric; otherwise, select the dual grid dominance relation, cosine similarity metric, and Pareto optimal subspace learning strategy. Obtain the non-dominated level of the current solution set by non-dominated sorting based on the corresponding dominance relation and store it in the convergence metric CM. Store the metric value of the current solution set in the diversity metric DM.

[0043] Step 6.2: Select candidate solutions for the next generation population based on the convergence metric CM obtained in Step 6.1. Add the solutions to the candidate solution set of the next generation population in ascending order of CM value until saturation or overflow occurs.

[0044] Step 6.3: Match the applicable maintenance strategy according to the evolutionary stage obtained in Step 3. If the evolution is in the first stage, truncate the candidate solution set of the next generation population in descending order according to the DM value obtained in Step 6.1. The retained part becomes the next generation population as the result of natural selection. At the same time, the non-dominated solutions in the next generation population are stored in the cumulative archive as input for the Pareto optimal subspace learning strategy. Otherwise, assign a solution with the maximum cosine similarity value to each reference vector in turn. The assigned solution becomes the next generation population as the result of natural selection.

[0045] The beneficial effects of this invention are: (1) The time-window constrained vehicle path planning method based on multi-objective evolutionary algorithm of this invention is a balance separation control mechanism based on a two-stage strategy. For the search behavior of each stage, corresponding components are used to guide the search direction, weaken the competition between convergence and diversity, and thus resolve the conflict characteristics between objectives to prevent the transportation scheme obtained under the background of complex problems from deteriorating severely on a certain objective, resulting in low utility performance. (2) The time-window constrained vehicle path planning method based on multi-objective evolutionary algorithm of this invention establishes a Kondratiev wave two-stage strategy model to carry the balance separation control mechanism. The Kondratiev wave two-stage strategy model has two different search tendencies and relatively regular dynamic periodic characteristics. The two originally independent stages will be executed alternately through automatic configuration components. This model can promote the information inheritance between stages and weaken the adverse effects of changes in a single state in each stage on the balance. (3) The time-window constrained vehicle path planning method based on multi-objective evolutionary algorithm of this invention proposes two enhanced domination methods that match the stage characteristics. The fusion of decomposition methods and Pareto advantages can ensure the effectiveness of distance-type and angle-type metrics in complex contexts, which is beneficial to balanced separation control. The first-stage hybrid sub-region domination method uses a boundary cross method with small parameter penalties to guide the convergence of the solution set. At the same time, it uses a distribution correction strategy with normalization as the main method and standardization as the auxiliary method to maintain an accurate distribution state and ensure that diversity is not deteriorated while emphasizing convergence, so that the solution can effectively approximate the real Pareto front. The second-stage dual grid domination method eliminates the cumbersome comparison mechanism within the grid and limits the occurrence of multiple pairs of equivalent solutions within the superbox by adaptively adjusting the superbox size, thereby providing dynamic selection pressure and reducing computational complexity. The Level metric used in the dual grid domination relationship in the second-level grid system plays a complementary role to the GS metric, implying a relaxed diversity representation while ensuring convergence. (4) The time window constrained vehicle path planning method based on multi-objective evolutionary algorithm of the present invention designs two environmental selection components that match the enhanced domination relationship to maintain the diversity of the first stage and explicitly improve the diversity of the second stage. In the first stage, which prioritizes convergence, the environment selection component employs a traditional truncation method to prevent excessive screening from severely degrading diversity. In the second stage, which prioritizes diversity, the environment selection component uses an allocation method to allow solutions that would otherwise be eliminated to enter the next generation. Since the reference vector set and the solution set have the same size, the closer the allocated solution is to the reference vector, the more uniformly it maintains a similar distribution. Simultaneously, the Pareto optimal subspace learning strategy is implemented before non-dominated sorting to find sub-regions that approximate the true Pareto optimal subspace and reduce the regions outside of it to compensate for the convergence loss. Attached Figure Description

[0046] Figure 1This is a logic block diagram of a time-window-constrained vehicle path planning method based on a multi-objective evolutionary algorithm according to the present invention.

[0047] Figure 2 This is a flowchart of the TS / KW-MaOEA used in the time window-constrained vehicle path planning method based on a multi-objective evolutionary algorithm of the present invention.

[0048] Figure 3 This is a graph showing the change in the average number of violating customers in the solution set obtained by TS / KW-MaOEA and other optimization algorithms on the Solomon-R101 test case.

[0049] Figure 4 This is a graph showing the change in the average transport mileage of the solution set obtained by TS / KW-MaOEA and other optimization algorithms on the Solomon-R101 test case.

[0050] Figure 5 This is a graph showing the change in the average number of vehicles dispatched from the solution set obtained by TS / KW-MaOEA and other optimization algorithms on the Solomon-R101 test case.

[0051] Figure 6 This is a graph showing the change in the average comprehensive cost of the solution set obtained by TS / KW-MaOEA and other optimization algorithms on the Solomon-R101 test case.

[0052] Figure 7 It is a transportation route map represented by a random solution in the set of optimal solutions obtained by TS / KW-MaOEA on the Solomon-R101 test case. Detailed Implementation

[0053] The present invention will now be described in detail with reference to the accompanying drawings and specific embodiments.

[0054] This invention discloses a time-window-constrained vehicle path planning method based on a multi-objective evolutionary algorithm, the process of which is as follows: Figure 1 As shown, the specific steps are as follows:

[0055] Step 1, Population Initialization

[0056] An initial population is obtained by random traversal insertion initialization method. The initial population represents the initial solution set of the time window-constrained vehicle path planning problem. The chromosome of any solution in the initial solution set represents the decision vector of the solution.

[0057] The population contains N solutions with chromosome length D, where D = cus_num + max_veh-1;

[0058] Where cus_num is the number of customers and max_veh is the maximum number of vehicles that can be dispatched;

[0059] The chromosome encoding of any solution x in the population is a randomly generated path planning scheme, x = {x1, x2, ..., x...} D};

[0060] Where, x i ...x j This represents the customer sequence of a transportation route service; specifically, x i =0 indicates a vehicle service center, x i >1+cus_num indicates the vehicle number; different transport routes are separated by the dispatched vehicle number;

[0061] The chromosome of any solution x can be decoded into configuration information of a transportation plan, including the number of dispatched vehicles, vehicle routes, vehicle distances, customer service sequences, the number of customers whose service time windows are violated, and the number of routes that violate customer needs; the chromosome of any solution x forms the genetic phenotype f(x) of that solution through function mapping, f(x) = {f1(x), f2(x), ..., f...} M (x)};

[0062] Among them, f i (x) represents the value obtained by transportation scheme x on the i-th objective function, and M represents the number of objective functions;

[0063] Step 2, Generate a set of reference vectors

[0064] The hyperplane construction method proposed by Das and Dennis is used to generate a set of reference vectors W, which contains N reference vectors to assist in the calculation of state characterization values ​​and the implementation of control mechanisms.

[0065] Step 3, Automatic Configuration of the Kondratiev Wave Two-Stage Strategy Model

[0066] The initial population undergoes iterative evolution. The evolutionary stage of the current population is obtained through the established Kondratiev wave two-stage strategy model, and the search tendency is adaptively adjusted. The timing of the search tendency switching is determined by the configuration to ensure the effectiveness of the balanced separation.

[0067] Specifically:

[0068] Obtain the number of objective functions and the current generation, and calculate the current evolutionary stage using formula (1):

[0069]

[0070] Where I represents the current evolutionary stage, G doneG represents the current generation number. max M represents the maximum number of generations, and M represents the number of objective functions. max With M min represents the maximum and minimum values ​​of the objective function in the statistical experiment, respectively, and rand(1) represents a random number in the interval [0, 1].

[0071] Step 4, Stage-based mating selection operation

[0072] Based on the evolutionary stage obtained in step 3, the dominance relationship and characterization value that match the stage characteristics are selected to indicate the state of the solution set. A mating pool is generated through a binary tournament with convergence measure as the primary selection criterion and diversity measure as the secondary selection criterion.

[0073] Specifically:

[0074] Step 4.1: Match the applicable dominance relation and representation value according to the evolutionary stage obtained in Step 3. If the evolution is in the first stage, select the mixed sub-region dominance relation and d2 distance metric; otherwise, select the dual grid dominance relation and SD density metric. Obtain the non-dominance level of the current solution set through the non-dominance sorting method based on the corresponding dominance relation and store it in the convergence metric CM. Store the metric value of the current solution set in the diversity metric DM.

[0075] The d1 and d2 distances used in the dominance relationship of the mixed subregions are calculated using formulas (4) and (5):

[0076]

[0077]

[0078] Where f(x) is a vector composed of objective function values, w represents any reference vector, and d1(x) and d2(x) represent the projected distance and perpendicular distance of the solution x to the reference vector w, respectively.

[0079] In the dual-grid dominance relationship, the construction of the first-layer grid system is based on formulas (6) and (7):

[0080] lb k =fmin k -(fmax k -fmin k ) / (2×div) (6)

[0081] ub k =fmax k +(fmax k -fmin k ) / (2×div) (7)

[0082] Among them, ub k With lb k fmax represents the upper and lower bounds of the grid in the k-th target, respectively. k with fmin k represents the maximum and minimum values ​​of the solution set in the k-th target, respectively, and div represents the number of segments in each target dimension;

[0083] In the dual-mesh dominance relationship, the hyperbox coordinates are assigned using formulas (8) and (9):

[0084] d k =(ub k -lb k ) / div (8)

[0085]

[0086] Where, d k G represents the width of the superbox in the k-th target. k (x) represents the coordinates of solution x in the k-th objective, f k (x) represents the k-th objective function value of the solution x;

[0087] Step 4.2: Randomly select solution a and solution b from the current solution set and compare them according to the convergence metric CM and diversity metric DM obtained in Step 4.1. The solution with the smaller CM value wins the comparison; if the CM values ​​are the same, the solution with the smaller DM value wins the comparison; otherwise, randomly select one of the two to win. Add the winning solution to the mating pool and repeat the above comparison operation until the mating pool is saturated.

[0088] Step 5: By applying ordered crossover and inversion mutation operators to the mating pool obtained in Step 3, a candidate offspring set is generated;

[0089] The solutions p1 and p2 obtained from the mating pool are used as the parent, with P c With a crossover probability of 0.9, several alleles in the father's chromosome are randomly selected as the start and end points for cutting. The chromosome fragment cut from father p1 is retained as genetic information in the allele sequence of the offspring. The remaining chromosome fragment from father p2 is added sequentially to the missing gene loci in the offspring. m With a mutation probability of 0.05, a gene sequence is randomly selected from the chromosomes of the offspring and its order is reversed.

[0090] Step 6, Stage-based environment selection operation

[0091] Based on the evolutionary stage obtained in step 3, select the dominance relationship and representation value that match the stage characteristics to indicate the state of the solution set, and select the next generation population from the union of the current population and the candidate offspring set obtained in step 5 using a matching maintenance strategy.

[0092] Specifically:

[0093] Step 6.1: Match the applicable dominance relation, representation value, and maintenance strategy according to the evolutionary stage obtained in Step 3. If the evolution is in the first stage, select the mixed sub-region dominance relation and d2 distance metric; otherwise, select the dual grid dominance relation, cosine similarity metric, and Pareto optimal subspace learning strategy. Obtain the non-dominated level of the current solution set by non-dominated sorting based on the corresponding dominance relation and store it in the convergence metric CM. Store the metric value of the current solution set in the diversity metric DM.

[0094] Cosine similarity is calculated using formula (10):

[0095]

[0096] Where similarity represents the cosine similarity value, and w represents any reference vector;

[0097] Step 6.2: Select candidate solutions for the next generation population based on the convergence metric CM obtained in Step 6.1. Add the solutions to the candidate solution set of the next generation population in ascending order of CM value until saturation or overflow occurs.

[0098] Step 6.3: Match the applicable maintenance strategy according to the evolutionary stage obtained in Step 3. If the evolution is in the first stage, truncate the candidate solution set of the next generation population in descending order according to the DM value obtained in Step 6.1. The retained part becomes the next generation population as the result of natural selection. At the same time, the non-dominated solutions in the next generation population are stored in the cumulative archive as input for the Pareto optimal subspace learning strategy. Otherwise, assign a solution with the maximum cosine similarity value to each reference vector in turn. The assigned solution becomes the next generation population as the result of natural selection.

[0099] The Pareto optimal subspace learning strategy is based on principal component analysis (PCA), using archive A to construct the PCA sample matrix K, where K∈R. |A|×D Each row of matrix K represents a solution x, and each column of decision variables represents a component of the solution set. Calculate the eigenvalues ​​of the matrix union K′ obtained by subtracting the mean from matrix K. Select decision variables with contribution coefficients ∈ > 0.99 as principal components. Retain the principal components and use the median of the non-principal components to cover the upper and lower bounds of their respective decision variables. The space formed by the new upper and lower bounds of all decision variables is the Pareto optimal subspace.

[0100] Step 7: Determine if the termination condition is met, that is, determine if the number of generations of the next generation obtained in Step 6 has reached the set maximum number of generations. If the maximum number of generations has not been reached, it means that the termination condition is not met, and jump to Step 3. Otherwise, output the next generation obtained in Step 6 as the optimal solution set.

[0101] The decision-maker selects a solution from the set of optimal solutions as the final decision vector, and obtains the transportation scheme configuration information corresponding to that vector by decoding.

[0102] The effects of this invention are further illustrated by the following simulation experiments.

[0103] First, given the test function, this experiment uses the R series (Solomon-R101) of the Solomon test set (Solomon test set) specifically for time-window-constrained vehicle routing planning problems. The Solomon-R101 data includes 100 customers and their two-dimensional location coordinates on a vehicle service center map, 100 sets of service time windows corresponding to each customer, one set of service center operating hours, and 100 loading and unloading requests corresponding to each customer. The customer location coordinates in Solomon-R101 simulate a widely random distribution under complex backgrounds, with no obvious clustering.

[0104] Then, a generalized mathematical model for the time-window-constrained vehicle path planning problem is given.

[0105] The transportation route map is defined as G(V, A), where V = {0, 1, ..., n} represents all nodes in the graph, and node 0 represents a vehicle service center; A is the set of all path arcs in the transportation route map. The unit transportation cost for path arc (i, j) is c. ij The transportation time is t ij Let the customer set be c = {1, ..., n}, and the customer's loading and unloading demand be q. i The customer's requested service time window is [e i , l i The vehicle arrives at the customer's location at time a. i The service duration is serv i Let the set of vehicles be K, and the maximum load capacity of the vehicles be... variable Indicate whether the k-th car serves path arc (i, j). If so, otherwise

[0106] The time window constraint in vehicle routing problems is divided into two types: hard time windows and soft time windows. Hard time windows require the vehicle to provide service within the customer's requested time window; otherwise, the customer will refuse the service. Soft time windows allow for different costs. Providing service earlier than the requested time incurs additional service costs for employees and negotiation costs with the customer, while delaying service incurs additional compensation costs for the customer. The biggest difference between hard and soft time windows lies in using penalty costs instead of waiting or refusal. and These represent the penalty coefficients for early service and delayed service, respectively.

[0107] In this invention, the optimization objectives of the time-window-constrained vehicle routing problem are threefold: minimizing the number of dispatched vehicles, minimizing the total distance of the transportation route, and calculating the comprehensive transportation cost using formula (11).

[0108]

[0109] The constraint of formula (12) guarantees that all customers can receive service:

[0110]

[0111] The constraint of formula (13) guarantees that each customer can only be served once:

[0112]

[0113] The constraint of formula (14) ensures that all transport routes originate from the vehicle service center:

[0114]

[0115] The constraint of formula (15) ensures that all transport routes end at the vehicle service center:

[0116]

[0117] Formula (16) guarantees the order in which adjacent nodes begin serving:

[0118]

[0119] Formula (17) represents the vehicle's load capacity constraint:

[0120]

[0121] Formula (18) represents the vehicle's load capacity constraint:

[0122]

[0123] Finally, the algorithm parameters are initialized, with the unit of measurement simplified to 1. The maximum number of dispatchable vehicles is set to 25, and the maximum load capacity of each vehicle is set to 200. One hundred 124-dimensional solutions are uniformly and randomly initialized in the search space as the initial population, i.e., N=100, D=124. The maximum number of generations is set to 100, the simulation is run 10 times, and the optimal solution set is output.

[0124] Simulation results are shown below Figures 3-7 This document records the evolutionary populations of this invention and three other commonly used algorithms for solving time-window-constrained vehicle routing problems (NSGA-III, NSGA-II-RPD, and IBEA). It decodes all solutions in the populations and calculates relevant information, such as... Figures 1-4 The following conclusions can be drawn:

[0125] This invention achieves significant advantages in minimizing overall transportation costs and the total mileage of transportation routes, and yields good results in minimizing the number of vehicles dispatched. Compared to the algorithm used in this invention, other algorithms achieve a slight advantage in minimizing the number of vehicles dispatched, but severely deteriorate the decision-making results of the transportation scheme on other objectives, even worsening the minimization of the average overall transportation cost by several times. Therefore, it can be concluded that the transportation scheme obtained through this invention has excellent utility performance.

[0126] A solution is randomly selected from the optimal solution combination, and the transportation scheme configuration information represented by this solution is obtained through decoding as follows:

[0127] Total transportation cost: 677,438.7483; Number of vehicles used: 17; Total distance of transportation routes: 3,306.0827.

[0128] Transportation route 1, 0 → 54 → 59 → 0;

[0129] Transportation route 2,0 → 44 → 40 → 61 → 51 → 53 → 41 → 0;

[0130] Transportation route 3,0 → 56 → 55 → 34 → 37 → 32 → 35 → 0;

[0131] Transportation route 4,0 → 39 → 47 → 48 → 0;

[0132] Transportation route 5,0 → 52 → 17 → 13 → 25 → 0;

[0133] Transportation route 6,0 → 6 → 3 → 10 → 4 → 1 → 100 → 0;

[0134] Transportation route 7,0 → 45 → 49 → 46 → 50 → 43 → 57 → 60 → 20 → 0;

[0135] Transportation route 8,0 → 27 → 33 → 28 → 36 → 31 → 29 → 42 → 0;

[0136] Transportation route 9,0 → 62 → 64 → 63 → 2 → 5 → 12 → 11 → 7 → 0;

[0137] Transportation route 10, 0 → 14 → 15 → twenty one → twenty three → 16 → 30 → 38 → 58 → 0;

[0138] Transportation route 11, 0→ 69 → 65 → 71 → 67 → 66 → 68 → 70 → 0;

[0139] Transportation route 12,0 → 73 → 85 → 79 → 81 → 78 → 84 → 74 → 75 → 76 → 77 → 80 → 82 → 83 → 0;

[0140] Transportation route 13,0 → 87 → 86 → 88 → 90 → 94 → 0;

[0141] Transportation route 14,0 → 96 → 91 → 89 → 92 → 95 → 98 → 93 → 0;

[0142] Transportation route 15, 0 → 99 → 8 → 9 → 0;

[0143] Transportation route 16,0 → 97 → 0;

[0144] Transportation route 17,0 → 72 → 19 → 18 → twenty two → 26 → twenty four → 0. The transportation route is shown in the attached diagram. Figure 5 .

Claims

1. A time window constrained vehicle routing method based on multi-objective evolutionary algorithm, characterized in that, The method is implemented according to the following steps: Step 1, an initial population is obtained by using a random traversal insertion initialization method, the initial population represents an initial solution set of the time window constrained vehicle routing problem, and a chromosome of any solution in the initial solution set represents a decision vector of the solution; Step 2, a reference vector set is generated by using a hyperplane construction method, the reference vector set contains the same number of reference vectors as the number of solutions in the initial population, and is used to assist in the calculation of state representation values and the implementation of a control mechanism; Step 3, the initial population is subjected to cyclic iterative evolution, an established two-stage strategy model of Kanpo is used to obtain an evolution stage in which an evolution generation number of the current population is located, and a search tendency is adaptively adjusted, and a time for switching the search tendency is configured to ensure the effectiveness of the balance separation; Step 4, a domination relation and a representation value that are suitable for characteristics of the evolution stage obtained in step 3 are selected to indicate a state of a solution set, and a mating pool is generated through a binary tournament selection with a convergence measure as a primary selection standard and a diversity measure as a secondary selection standard; Step 5, a candidate offspring set is generated by using an ordered crossover and a reversal mutation operator on the mating pool obtained in step 4; Step 6, a domination relation and a representation value that are suitable for characteristics of the evolution stage obtained in step 3 are selected to indicate a state of a solution set, and a maintenance strategy that matches the state is selected to screen a next generation population from a union set of the current population and the candidate offspring set obtained in step 5; Step 7, whether the evolution generation number of the next generation population obtained in step 6 reaches a set maximum evolution generation number is judged, if the maximum evolution generation number is not reached, it is indicated that a termination condition is not met, and the method jumps to step 3, otherwise, the next generation population obtained in step 6 is output as an optimal solution set.

2. The multi-objective evolutionary algorithm based time window constrained vehicle routing method according to claim 1, characterized in that, The two-stage strategy model of Kanpo in step 3 includes a balance separation control mechanism, and contains a first stage, a second stage and a third stage; in the first stage, a balance state of a solution set is separated into a convergence state and a distribution state, and a corresponding compromise search tendency is separated into convergence emphasis and diversity emphasis, and stages with different search tendencies are alternately executed according to a switching rule of an automatic configuration component; in the second stage, the stages with different search tendencies correspondingly use different emphasis domination rules, selection standards and maintenance strategies in a domination relation, a mating selection component and an environment selection component; in the third stage, the convergence state and the distribution state of different stages are represented by using measurement values that match characteristics of the stages. 3.The multi-objective evolutionary algorithm based time window constrained vehicle routing method according to claim 2, characterized in that, The automatic configuration component of the two-stage strategy model of Kanpo in step 3 is based on statistical experiments; a switching rule is obtained by simulating a macroeconomic regulation and control method of the Keynesian school and analyzing a ratio parameter μ of a first stage evolution generation number and a maximum evolution generation number, in at least a first 75% of a search process, the automatic configuration component guides the search to stay in the first stage with convergence emphasis, and in at most a last 25% of the search process, the automatic configuration component guides the two stages to alternately execute according to a number of objective functions and a certain randomness.

4. The method of claim 3, wherein, In step 3, a number of objective functions and a current evolution generation number are obtained, and an evolution stage in which the current evolution generation number is located is calculated through formula (1): where I denotes the evolutionary stage at the current evolutionary generation, G done denotes the current evolutionary generation, G max denotes the maximum evolutionary generation, M denotes the number of objective functions, M max and M min denote the maximum and minimum values of the number of objective functions in the statistical experiment, respectively, and rand(1) denotes a random number in the interval [0, 1].

5. The multi-objective evolutionary algorithm based time window constrained vehicle routing method according to claim 4, characterized in that, The dominance relationship matching the phase characteristics in step 4 is a hybrid sub-region dominance relationship and a double-grid dominance relationship respectively: Regarding the mixed sub-region dominance relation: Assume that u and v are two solutions, u dominates v, denoted as 1) 2) u and v are equivalent in the Pareto sense; a) u and v are associated with the same reference vector after normalization; 1.

2. But under standardization, u and v are related to the same reference vector; PBI(u) < PBI(v), dl(u) < dl(v) and b) u and v are associated with different reference vectors when normalized, but the same reference vector when standardized; and PBI(u) < PBI(v), and dl(u) < dl(v), wherein the notation indicates the normalized metric value, the notation indicates Pareto dominance; di and d2 indicate the projected and perpendicular distances of the solution on the reference vector associated therewith; SD indicates the sub-region density of the solution, i.e. the number of solutions associated with the same reference vector as the solution; PBI indicates the metric value obtained using the penalty-based boundary crossover used in the MOEA / D algorithm and formula (2): PBI(x) = d1(x) + θd2(x), θ = 1 (2) Regarding the dual grid domination relation: Assume u and v are two solutions, u dominates v, denoted as 1) 2) u and v are equivalent in the Pareto sense; a) under a first level grid system, b) u and v are equivalent in the first-level grid system; GS(u) < GS(v) and Level(u) < Level(v); wherein The notation indicates the grid dominance used in the GrEA algorithm, i.e. the first level grid hierarchy; GS indicates the sum of the coordinate parameters of the hypercube to which the solution belongs, solutions within the same hypercube share the GS metric value; Level indicates the horizontal layer ordering level of the hypercube, hypercubes at the same horizontal layer share the Level metric value.

6. The multi-objective evolutionary algorithm based time window constrained vehicle routing method according to claim 5, characterized in that, The double-grid dominance relationship in step 4 realizes self-adaptive adjustment of the hyper-box size through formula (3): where (||·||) denotes the number of interior elements, div denotes the number of cuts for uniform division of each target axis, a is a valve factor related to the number of objectives M, a = 2M(M - 2), P t denotes the tth generation population, PF1 denotes the set of non-dominated solutions in the tth generation population on the first layer Pareto front.

7. The multi-objective evolutionary algorithm based time window constrained vehicle routing method according to claim 6, characterized in that, The step 4 is specifically: Step 4.1, according to the evolution phase obtained in step 3, match the applicable dominance relationship and the characteristic value, if the evolution is in the first phase, the hybrid sub-region dominance relationship and d2 distance measure are selected, otherwise the double-grid dominance relationship and SD density measure are selected; the non-dominated level of the current solution set is obtained through the non-dominated sorting method based on the corresponding dominance relationship and stored in the convergence measure CM, and the measure value of the current solution set is stored in the diversity measure DM; Step 4.2, randomly select solution a and solution b from the current solution set and compare them according to the convergence measure CM and the diversity measure DM obtained in step 4.1, wherein the solution with smaller CM value wins in the comparison; if the CM values of the two are the same, the solution with smaller DM value wins in the comparison; otherwise, one of them is randomly selected to win; the winning solution is put into the mating pool, and the above comparison operation is repeated until the mating pool is saturated.

8. The multi-objective evolutionary algorithm based time window constrained vehicle routing method according to claim 7, characterized in that, Step 6 is specifically: Step 6.1, according to the evolution phase obtained in step 3, match the applicable dominance relationship, characteristic value and maintenance strategy, if the evolution is in the first phase, the hybrid sub-region dominance relationship and d2 distance measure are selected, otherwise the double-grid dominance relationship, cosine similarity measure and Pareto optimal subspace learning strategy are selected; The non-dominated level of the current solution set is obtained through the non-dominated sorting method based on the corresponding dominance relationship and stored in the convergence measure CM, and the measure value of the current solution set is stored in the diversity measure DM; Step 6.2, according to the convergence measure CM obtained in step 6.1, select the candidate solutions of the next generation population, and put the solutions into the candidate solution set of the next generation population layer by layer according to the ascending order of CM value until saturation or overflow; Step 6.3, according to the evolution phase obtained in step 3, match the applicable maintenance strategy, if the evolution is in the first phase, the DM value obtained in step 6.1 is used to truncate the candidate solution set of the next generation population in descending order, the retained part becomes the next generation population as a result of natural selection, and the non-dominated solutions in the next generation population are stored in the cumulative archive as input for the Pareto optimal subspace learning strategy; Otherwise, assign each reference vector to a solution with the maximum cosine similarity value, and the assigned solution becomes the next generation population as a result of natural selection.

Images