CN114707693A - A method and system for vehicle distribution path planning - Google Patents

Abstract

The invention relates to a vehicle distribution path planning method and system. First, the node location information, customer information and vehicle information of the distribution center are acquired, and an objective function is established based on the shortest total distribution vehicle mileage; then the scanning algorithm is used to obtain An initial solution takes the objective function as the target value, and uses the improved large-scale neighborhood search algorithm combined with the Metropolis criterion to select the optimal vehicle delivery route. The use of the improved large-scale neighborhood search algorithm can filter out the optimal distribution path for vehicles, so as to formulate a more reasonable distribution path scheme, improve the distribution efficiency, and at the same time greatly reduce the total vehicle distribution mileage and reduce the manual formulation of the scheme. Work intensity.

Description

A method and system for vehicle distribution path planning

technical field

The invention relates to the field of goods distribution, in particular to a large-scale vehicle distribution path planning method and system based on an improved adaptive large-scale neighborhood search algorithm.

Background technique

The research of vehicle routing is a hot issue in the past ten years, and it has become a development trend of logistics distribution. It is closely related to resource distribution and path planning in real life, and mainly solves the problem of how to reasonably plan and determine the optimal vehicle path when multiple delivery vehicles with limited capacity serve multiple customer nodes. With the increase in the scale of the distribution problem, that is, the increase in the number of distribution customer nodes, the formulation of the distribution scheme becomes more complicated. The existing schemes generally rely on manual experience, but it is difficult to obtain a better vehicle path only by relying on manual experience. The distribution plan, resulting in a waste of labor and time.

Therefore, there is an urgent need in the art for a planning solution for the optimal distribution path of vehicles that can be used for large-scale customers.

SUMMARY OF THE INVENTION

The purpose of the present invention is to provide a vehicle distribution path planning method and system, select the distribution limit of the maximum vehicle load, take the vehicle routing problem with capacity constraints as the model of vehicle distribution, and take the shortest distribution mileage of each vehicle in the distribution scheme as To achieve the goal, an improved large-scale neighborhood search algorithm is introduced to select the optimal distribution route, so as to formulate a reasonable distribution route plan and reduce the work intensity of manual plan formulation.

For achieving the above object, the present invention provides the following scheme:

A vehicle distribution path planning method, the method comprising:

Acquiring the location information, customer information and vehicle information of the distribution center node; the customer information includes: the number of customers, the demand for goods of each customer, and the location of each customer node; the vehicle information includes: the number of vehicles and the maximum load information of the vehicle;

According to the distribution center node location information, customer information and vehicle maximum load information, establish an objective function with the shortest total distribution vehicle mileage as the goal; the constraint condition of the objective function is that each customer is served by one and only one vehicle and each The total customer demand served by the vehicle does not exceed the vehicle load;

determining an initial coding sequence; the initial coding sequence is a coded representation of the initial delivery route, and the initial coding sequence includes each client node on the initial delivery route and the delivery sequence of each client node;

Calculate the target value of the initial coding sequence according to the target function;

Calculate the initial temperature of the Metropolis criterion according to the target value of the initial coding sequence, and start the optimal delivery route iteration;

The optimal delivery route iteration specifically includes:

According to the weight of each removal operator and each insertion operator, use the roulette method to select a removal operator and an insertion operator, and according to the selected removal operator and insertion operator, respectively, for the current solution Perform an operator removal operation and an operator insertion operation to obtain a local solution; the removal operators include: random customer node removal, random subpath removal, similar customer node removal and ring similarity removal; all The insertion operator includes: cost greedy insertion and regret value insertion; the current solution in the first iteration process is the initial coding sequence;

If the target value of the local solution is smaller than the target value of the current solution, the local solution is taken as the current solution; if the target value of the local solution is greater than the target value of the current solution, then according to the Metropolis criterion, combined with the temperature value of this iteration, the probability Accept the local solution and update the temperature value; the temperature value in the first iteration process is the initial temperature value;

According to the usage times and preset scores of each insertion operator and each removal operator, the corresponding operator weights are updated, and the update interval is the preset number of iterations;

Until the target value of the local solution is not less than the current solution iteration times and reaches the preset number of times, the iteration is stopped, and the vehicle optimal coding sequence is output.

其中，F表示目标函数；m为配送能够使用的最大车辆数；n为客户数量；d_ij表示任意两节点i和j之间的距离；x_ijk为0-1规划，表示车辆k是否经过节点i和j，0表示车辆k不经过节点i和j，1表示车辆k经过节点i和j。In some embodiments, the objective function aiming at the shortest total delivery vehicle mileage is:

Among them, F represents the objective function; m is the maximum number of vehicles that can be used for delivery; n is the number of customers; _{d ij} represents the distance between any two nodes i and j; i and j, 0 means that vehicle k does not pass through nodes i and j, and 1 means that vehicle k passes through nodes i and j.

In some embodiments, the determining the initial coding sequence specifically includes:

Taking the position of the distribution center node as the pole, and taking the ray of the distribution center node and the customer node closest to the distribution center node as the polar axis, a polar coordinate system is established, and the position coordinates of the distribution center node and the customer node are converted. is polar coordinates;

Sort client nodes in order of polar angle size;

On the premise that the constraints are met, each customer node is inserted into the solution in turn according to the objective function; the solution includes several distribution paths; during the insertion process, if a customer node that does not meet the constraints is encountered, and all If the client node that does not meet the constraints is not the last client node, the client node that does not meet the constraints is skipped;

Search whether there is an uninserted customer node j, and if so, calculate the insertion cost c _ij between each two customer nodes inserted by the customer node j into each delivery route, and insert the two corresponding to the lowest cost. between client nodes; where c _ij =d _i,j +d _j,i+1 -d _i,i+1 ; in the formula, d _i,j represents the distance between client node i and client node j, d _j,i+1 represents the distance between client node j and client node i+1, and d _i,i+1 represents the distance between client node i and client node i+1;

Traverse all uninserted client nodes until all client nodes are inserted into the solution to obtain the initial encoding sequence.

In some embodiments, calculating the initial temperature of the Metropolis criterion according to the target value in the initial coding sequence specifically includes:

所述初始温度T₀表示在劣解接受概率为η的情况下，要接受比当前目标值f(x)差β％的劣解，所对应的温度为T₀。Calculate the initial temperature T ₀ of the Metropolis criterion according to the target value in the initial coding sequence;

The initial temperature T ₀ indicates that under the condition that the acceptance probability of the inferior solution is η, the inferior solution that is worse than the current target value f(x) by β% is to be accepted, and the corresponding temperature is T ₀ .

In some embodiments,

The random customer node removal specifically includes: if there is an infeasible sub-path solution in the current solution, randomly removing the customer nodes in the infeasible sub-path until the sum of the customer demand does not exceed the maximum vehicle load; If the number of customer nodes does not exceed the removal upper limit, select random customer nodes for removal until the number of removed customer nodes reaches the removal upper limit; the solution of the infeasible sub-path is that the sum of the customer demand of the sub-path is greater than the maximum vehicle load; the removal upper limit is generated by a random function;

The random sub-path removal specifically includes: if there is an infeasible sub-path solution in the current solution, removing the infeasible sub-path; if the number of removed client nodes does not exceed the removal upper limit, selecting a random sub-path for removal. Remove until the number of removed client nodes reaches the removal limit;

The removal of similar client nodes specifically includes: randomly selecting a client node i in the current solution, calculating the similarity between the client node i and other client nodes j, and sequentially removing the client node with the largest similarity until the client node is removed. The number of nodes reaches the removal limit;

The removal of the similarity of the ring area specifically includes: calculating the similarity based on three spatial relationships between the ring areas where each sub-path in the current solution is located;

The similarity between the two sub-paths is:

Among them, σ ₁ >σ ₂ >σ ₃ , s _i and _{s j} represent the annular area corresponding to the two sub-paths ri and r _j in the current solution, then the similarity coefficient of a single sub-path _ri is:

中选出相似度系数最大的子路径r_m镠x，通过轮盘赌的方法选出非最大相似度系数的子路径r_d，之后移除r_m镠x和r_d中的所有客户节点，

表示所有子路径的相似度系数集合，l表示当前解中子路径总数；exist

Select the sub-path r _{m x x} with the largest similarity coefficient, select the sub-path r _d with the non-maximum similarity coefficient through the roulette method, and then remove all customer nodes in r _{m x} and r _d ,

Represents the set of similarity coefficients of all subpaths, and l represents the total number of subpaths in the current solution;

The cost-greedy insertion specifically includes: inserting each removed client node into the sub-path of the current solution, and selecting the sub-path with the smallest insertion cost and satisfying the load restriction for each insertion;

The regret value insertion specifically includes: inserting each removed client node into the sub-path of the current solution, and selecting the sub-path with the largest regret value for each insertion; the regret value is: Δc=c _m2 -c _m1 , where c _m1 represents the smallest insertion cost and c _m2 represents the next smallest insertion cost.

In some embodiments, the update formula of the operator weight is:

表示算子q在τ次迭代内的使用次数，对应的预设评分为

在第y个迭代间隔内，第q个算子对应的算子权重为

在第y+1个迭代间隔内，第q个算子对应的算子权重为

μ∈[0,1]是控制算子影响的因子，该值越小，则第y个迭代间隔的权重对第y+1个迭代间隔的权重更新影响越大；若μ＝0，则权重完全取决于上一迭代间隔的权重；若μ＝1，则只考虑迭代间隔内评分，若0<μ<1，则同时考虑迭代间隔内的得分与之前的权重值。where τ represents the interval τ iterations, y represents the y-th iteration interval, q represents the q-th operator,

represents the number of times the operator q is used in τ iterations, and the corresponding preset score is

In the yth iteration interval, the operator weight corresponding to the qth operator is

In the y+1th iteration interval, the operator weight corresponding to the qth operator is

μ∈[0,1] is the factor that controls the influence of the operator. The smaller the value, the greater the influence of the weight of the y-th iteration interval on the weight update of the y+1-th iteration interval; if μ=0, the weight It depends entirely on the weight of the previous iteration interval; if μ=1, only the score in the iteration interval is considered, and if 0<μ<1, the score in the iteration interval and the previous weight value are considered at the same time.

In some embodiments, before the determining the initial coding sequence, the method further includes:

Chromosome coding is performed on the distribution center node and each of the customer nodes.

In some embodiments, the update formula of the temperature value is: T′=T*α

Among them, T′ represents the temperature after updating, T represents the current temperature before updating, α∈(0,1) is the cooling speed, and the value of α is a preset constant parameter.

In some embodiments, according to the Metropolis criterion, combined with the temperature value of this iteration, the probability of accepting a local solution includes:

其中，P_Δ的值与[0,1]间的随机数r，若P_Δ≥r，则接受劣解作为下一次迭代的当前解，T表示当前温度，f(x)表示当前解x的目标值，f(x′)表示局部解x′的目标值。According to the Metropolis criterion, combined with the temperature value of this iteration, the probability of accepting the local solution is

Among them, the value of P _Δ and the random number r between [0,1], if P _Δ ≥ r, the inferior solution is accepted as the current solution for the next iteration, T represents the current temperature, and f(x) represents the current solution x target value, f(x') represents the target value of the local solution x'.

The present invention also provides a vehicle distribution path planning system, the system includes:

The information acquisition module is used to acquire the node location information, customer information and vehicle information of the distribution center; the customer information includes: the number of customers, the demand for goods of each customer, and the location of each customer node; the vehicle information includes: the number of vehicles and Vehicle maximum load information;

The objective function determination module is used to establish an objective function with the shortest total distribution vehicle mileage as the objective according to the distribution center node location information, customer information and vehicle maximum load information; the constraint condition of the objective function is that each customer consists of and only Served by one vehicle and the total customer demand served by each vehicle does not exceed the vehicle load;

The initial solution building module is used to determine the initial coding sequence; the initial coding sequence is the coded representation of the initial distribution path, and the initial coding sequence includes each customer node on the initial distribution path and the distribution of each customer node order;

an initial solution target value calculation module, for calculating the target value of the initial coding sequence according to the target function;

an optimal delivery route iteration module, configured to calculate the initial temperature of the Metropolis criterion according to the target value of the initial coding sequence, and start the optimal delivery route iteration;

The optimal delivery route iteration specifically includes:

According to the weight of each removal operator and each insertion operator, use the roulette method to select a removal operator and an insertion operator, and according to the selected removal operator and insertion operator, respectively, for the current solution Perform an operator removal operation and an operator insertion operation to obtain a local solution; the removal operators include: random customer node removal, random subpath removal, similar customer node removal and ring similarity removal; all The insertion operator includes: cost greedy insertion and regret value insertion; the current solution in the first iteration process is the initial coding sequence;

If the target value of the local solution is smaller than the target value of the current solution, the local solution is taken as the current solution; if the target value of the local solution is greater than the target value of the current solution, then according to the Metropolis criterion, combined with the temperature value of this iteration, the probability Accept the local solution and update the temperature value; the temperature value in the first iteration process is the initial temperature value;

According to the usage times and preset scores of each insertion operator and each removal operator, the corresponding operator weights are updated, and the update interval is the preset number of iterations;

Until the target value of the local solution is not less than the current solution iteration times and reaches the preset number of times, the iteration is stopped, and the vehicle optimal coding sequence is output.

According to the specific embodiments provided by the present invention, the present invention discloses the following technical effects:

The invention provides a vehicle distribution path planning method and system. First, the node location information, customer information and vehicle information of the distribution center are acquired, and an objective function is established based on the shortest total distribution vehicle mileage; An initial solution is obtained, the target value of the initial coding sequence is determined by the objective function, and the optimal vehicle distribution route is selected by the improved large-scale neighborhood search algorithm combined with the Metropolis criterion. The present invention utilizes the improved large-scale neighborhood search algorithm to screen out the optimal distribution path for vehicles, thereby formulating a more reasonable distribution path scheme, improving distribution efficiency, and at the same time greatly reducing the total vehicle distribution mileage and reducing manual formulation. Work intensity of the program.

Description of drawings

In order to more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the accompanying drawings required in the embodiments will be briefly introduced below. Obviously, the drawings in the following description are only some of the present invention. In the embodiments, for those of ordinary skill in the art, other drawings can also be obtained according to these drawings without any creative effort.

FIG. 1 is a flowchart of a method for planning a vehicle distribution path provided by Embodiment 1 of the present invention.

FIG. 2 is a schematic diagram of the location of a distribution center and each node according to Embodiment 1 of the present invention.

FIG. 3 is a schematic diagram of a spatial relationship between two sub-paths corresponding to an annular area according to Embodiment 1 of the present invention.

FIG. 4 is a schematic diagram of a process for four vehicles to complete all delivery tasks and a schematic diagram of a client node that each vehicle needs to serve, according to Embodiment 1 of the present invention.

FIG. 5 is a block diagram of a vehicle distribution path planning system according to Embodiment 2 of the present invention.

Detailed ways

The technical solutions in the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings in the embodiments of the present invention. Obviously, the described embodiments are only a part of the embodiments of the present invention, but not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those of ordinary skill in the art without creative efforts shall fall within the protection scope of the present invention.

The purpose of the present invention is to provide a vehicle distribution path planning method and system, select the distribution limit of the maximum vehicle load, take the vehicle routing problem with capacity constraints as the model of vehicle distribution, and take the shortest distribution mileage of each vehicle in the distribution scheme as To achieve the goal, an improved large-scale neighborhood search algorithm is introduced to select the optimal distribution route, so as to formulate a reasonable distribution route plan and reduce the work intensity of manual plan formulation. First, obtain the basic information for determining the vehicle path; the basic information includes the geographical location of the distribution center node and each customer node, the number of customers, the customer's demand for goods and the maximum load of the vehicle; then establish a vehicle distribution path with the shortest mileage as the goal. Then, the initial solution of the model is determined by the scanning algorithm. Finally, using the objective function as the target value, an improved large-scale neighborhood search algorithm is used to determine the optimal vehicle distribution path.

In order to make the above objects, features and advantages of the present invention more clearly understood, the present invention will be described in further detail below with reference to the accompanying drawings and specific embodiments.

Example 1:

As shown in FIG. 1 , this embodiment provides a vehicle distribution path planning method, the method includes:

S1. Obtain the location information, customer information and vehicle information of the distribution center node; the customer information includes: the number of customers, the demand for goods of each customer, and the location of each customer node; the vehicle information includes: the number of vehicles and the maximum load information of the vehicle . The location of the distribution center and each customer node is shown in Figure 2.

S2. Establish an objective function with the shortest total distribution vehicle mileage as the goal according to the distribution center node location information, customer information and vehicle maximum load information; the constraint condition of the objective function is that each customer is served by one and only one vehicle And the total customer demand served by each vehicle does not exceed the vehicle load.

其中，F表示目标函数；m为配送能够使用的最大车辆数；n为客户数量；d_ij表示任意两节点i和j之间的距离；x_ijk为0-1规划，表示车辆k是否经过节点i和j，1表示车辆k经过节点i和j，0表示车辆k不经过节点i和j，0表示节点i和节点j点间不连通，即车辆k不经过节点i和节点j之间的路径，也就不经过该路段的两个端点，即节点i和节点j。The objective function with the shortest total delivery vehicle mileage as the goal is:

Among them, F represents the objective function; m is the maximum number of vehicles that can be used for delivery; n is the number of customers; _{d ij} represents the distance between any two nodes i and j; i and j, 1 means that the vehicle k passes through the nodes i and j, 0 means that the vehicle k does not pass through the nodes i and j, 0 means that there is no connection between the node i and the node j, that is, the vehicle k does not pass through the node i and the node j. The path does not pass through the two endpoints of the road segment, namely node i and node j.

Considering the distribution of distribution between vehicles and customers, as well as the load of the distribution vehicle itself, will affect the choice of vehicle distribution path, so in order to obtain a better distribution path, it is necessary to set certain constraints on the objective function. That is, the constraint condition of the objective function is that each customer is served by one and only one vehicle, and the total customer demand served by each vehicle does not exceed the vehicle load.

Then, chromosome coding is performed on the distribution center node and each of the customer nodes. In order to facilitate the distinction between the distribution center and the customer node, 0 is used to represent the distribution center, that is, the departure point of the delivery vehicle, and each customer node is represented by a natural number.

S3. Determine an initial coding sequence; the initial coding sequence is a coded representation of the initial delivery route, and the initial coding sequence includes each client node on the initial delivery route and the delivery sequence of each client node. The delivery sequence is the sub-path and sequence of delivery. For example, the initial coding sequence includes distribution centers (represented by 0 here) and ordered customer points (represented by natural numbers such as 1, 2, and 3). For example, for a vehicle routing problem with one distribution center and six distribution points, an initial coding sequence can be expressed as {0, 3, 4, 2, 0, 5, 1, 6, 0}, and {} is a solution The code means that the vehicle starts from the distribution center (0), passes through customer points 3, 4, and 2 in turn, and then returns to the distribution center, and then starts from the distribution center (0) and passes through customers 5, 1, and 6 and then returns to the distribution center. For distribution center 0, customer points 1, 2, 3, 4, 5, 6, after step S3, can form {0, 3, 4, 2, 0, 5, 1, 6, 0}.

The steps of determining the initial coding sequence include:

S31. Taking the position of the distribution center node as the pole, and taking the ray where the distribution center node and the customer node closest to the distribution center node are located as the polar axis, a polar coordinate system is established, and the positions of the distribution center node and the customer node are set as polar axes. Coordinates are converted to polar coordinates.

S32: Sort the client nodes according to the order of polar angle size. The ordering rule is from small to large or from large to small.

S33. On the premise of satisfying the constraints, insert each customer node into the solution in sequence according to the objective function; the solution includes several delivery paths. That is, after solving the distribution of a customer node, put the solved customer node and path into the solution.

During the insertion process, if a client node that does not meet the constraint condition is encountered, and the client node that does not meet the constraint condition is not the last client node, the client node that does not meet the constraint condition is skipped.

S34. Search whether there is an uninserted customer node j, and if so, calculate the insertion cost c _ij between each two customer nodes inserted by the customer node j into each delivery route in turn, and insert the lowest cost corresponding to the _between the _{two client nodes of} Distance, d _{j, i+1} represents the distance between client node j and client node i+1, and d _{i, i+1} represents the distance between client node i and client node i+1.

For example, {0, 3, 4, 2, 0, 5, 1, 6, 0}, the solution is just {0} at the beginning, and then consider the insertion of the points 1, 2, 3, 4, 5, 6 Cost, based on the comparison of the insertion cost to determine which point is placed after 0, such as inserting 3, it becomes {0, 3}, and so on, until {0, 3, 4, 2,}, and then inserting 5 when the point is found The sum of the customer requirements of 3, 4, 2, and 5 exceeds the load, so 5 cannot be inserted, so 0 is inserted first to become {0, 3, 4, 2, 0}, and then 5 is inserted to form {0, 3, 4, 2, 0, 5}, until all points are inserted, forming {0, 3, 4, 2, 0, 5, 1, 6, 0}.

S35: Traverse all uninserted client nodes until all client nodes are inserted into the solution, and obtain the initial coding sequence.

S4. Calculate the target value of the initial coding sequence according to the target function.

The target value refers to the total vehicle path length. Taking the solution {0, 3, 4, 2, 0, 5, 1, 6, 0} as an example, the target value f(x)=d ₀₃ +d ₃₄ +d ₄₂ +d ₂₀ +d ₀₅ +d ₅₁ + d ₁₆ +d ₆₀ , where d _ij is the distance between point i and point j. The target value here can also be understood as the fitness of a chromosome in the genetic algorithm, and its essence is to evaluate the quality of a solution.

S5. Calculate the initial temperature of the Metropolis criterion according to the target value of the initial coding sequence, and start the iteration of the optimal delivery route, and output the optimal coding sequence of the vehicle after the iteration.

所述初始温度T₀表示在劣解接受概率为η的情况下，要接受比当前目标值f(x)差β％的劣解，所对应的温度为T₀，其中η和β一般是人为规定，无固定范围，本实施例中令η＝50和β＝20。The initial temperature T ₀ of the Metropolis criterion is calculated as:

The initial temperature T ₀ indicates that under the condition that the acceptance probability of the inferior solution is η, the inferior solution that is worse than the current target value f(x) by β% is to be accepted, and the corresponding temperature is T ₀ , where η and β are generally artificial. It is specified that there is no fixed range, and η=50 and β=20 in this embodiment.

The optimal delivery route iteration specifically includes:

A1. According to the weight of each removal operator and each insertion operator, use the roulette method to select a removal operator and an insertion operator, and according to the selected removal operator and insertion operator, respectively The current solution is subjected to operator removal and operator insertion to obtain a local solution; the removal operators include: random client node removal, random subpath removal, similar client node removal, and loop similarity removal ; The insertion operator includes: cost greedy insertion and regret value insertion; the current solution in the first iteration process is the initial coding sequence. At the beginning: the weight of each operator (four types of removal operators and two types of insertion operators in this embodiment) is 1 at the beginning of the algorithm, and the weight is initially predetermined.

In step S3, an initial solution of the vehicle routing problem has been constructed. This initial solution can actually be used as the final answer to the problem, but the quality of the initial solution is generally poor. Therefore, the subsequent steps are mainly to optimize the initial solution. This optimization process It belongs to the main body of the adaptive large-scale neighborhood search algorithm (ALNS for short). The ALNS algorithm itself does not include a description of the method for generating the initial solution, so other methods need to be used to construct the initial solution. Here, the scanning algorithm is selected, and the method can also be the saving algorithm, the adjacent insertion algorithm, etc. The optimization idea of the ALNS algorithm can be understood as the continuous destruction and reorganization of the solution, and the solution will transform to a better direction during this process.

Just as in the simulated annealing algorithm the current temperature is calculated at each iteration, removal and insertion are inherent operations and core features of the ALNS algorithm. The optimization process of the ALNS algorithm can be simply described as: complete solution (like the initial solution, including complete path information) → execute the removal operator (remove some points or some sub-paths in the complete solution, and the removed solution It can be called a non-complete solution) → execute the insertion operator (re-insert the process described by the removed point removal operator into the non-complete solution, and form a complete solution again) → complete solution, thus forming a closed loop, continuously Iterative optimization. One of the removal operator and the insertion operator forms an iteration, and only one removal operator and insertion operator are selected for each iteration.

The selection of the removal operator (select one of multiple removal operators) and the selection of the insertion operator (select one of multiple insertion operators) are both roulette, and the selections of the two do not interfere with each other. .

Roulette selection method is a common method in algorithms, which is commonly used in the selection of chromosomes in genetic algorithms, and the selection of operators here is similar.

其中，w_u与该公式等号前的P_u对应，指的是第u个算子的权重(此处是为了说明个体概率的计算方式，因而以第u个算子为例)，w_g与公式中的

对应，指w_g中g从1开始至h为止的w_g之和，其中的g作用域仅在求和符号之内，指代求和范围符号，而其中h是该段落中的h，即移除算子的总数。算子的累计概率为

其中，P_g为个体概率，

为累计概率；个体概率与累计概率均为轮盘赌这一通用方法中的原有值。与上一个问题中的g同理，Pg在公式

中，P为个体概率，g为其下标，表示求和符号中的序列符号，即Pg为第g个算子的个体概率(由前面的个体概率公式求出)，其理论含义为第g个算子在所有算子(从第1至第h个算子)中，被选取的概率。

表示前n个算子的个体概率之和，如

表示算子1与算子2的概率和。生成一个随机数r∈[0,1]，从n＝1开始，依次比较r与

当第一次出现

时，选择过程结束，第u个算子即为选出的算子。Taking 1 removal operator selected from h removal operators as an example, set the corresponding weights of the _h removal operators as w ₁ , w ₂ , ..., wh , then the uth operator corresponds to the individual probability

Among them, w _u corresponds to P _u before the equal sign of the formula, and refers to the weight of the u-th operator (here is to explain the calculation method of individual probability, so the u-th operator is taken as an example), w _g with the formula

Corresponding, refers to the sum of w _g from 1 to h in w _g , where g scope is only within the summation symbol, refers to the summation range symbol, and where h is the h in the paragraph, that is The total number of removal operators. The cumulative probability of the operator is

Among them, P _g is the individual probability,

is the cumulative probability; both the individual probability and the cumulative probability are the original values in the general method of roulette. Similar to g in the previous question, Pg is in the formula

Among them, P is the individual probability, and g is the subscript, indicating the sequence symbol in the summation symbol, that is, Pg is the individual probability of the gth operator (calculated by the previous individual probability formula), and its theoretical meaning is the gth The probability that an operator is selected among all operators (from the 1st to the hth operator).

represents the sum of the individual probabilities of the first n operators, such as

represents the probability sum of operator 1 and operator 2. Generate a random number r∈[0,1], starting from n=1, compare r with

when it first appeared

When , the selection process ends, and the u-th operator is the selected operator.

The removal operator is executed first for the current solution, and then the insertion operator is executed. The removal operator is used to remove some customer points from the current solution until the removed customer points reach the predetermined removal upper limit, thereby obtaining an incomplete solution; the insertion operator is used to remove the removed customer points according to The rules described by the insertion operator are re-inserted into the incomplete solution, so as to obtain a new local solution, so as to achieve the purpose of searching for a better solution.

The random customer node removal specifically includes: if there is an infeasible sub-path solution in the current solution, randomly removing the customer nodes in the infeasible sub-path until the sum of the customer demand does not exceed the maximum vehicle load; If the number of customer nodes does not exceed the removal upper limit, select random customer nodes for removal until the number of removed customer nodes reaches the removal upper limit; the solution of the infeasible sub-path is that the sum of the customer demand of the sub-path is greater than the maximum vehicle load; The removal cap is generated by a random function. The upper limit of removal is generally a random number between 20% and 30% of the number of client nodes. The random number is generally generated by a random function in the program.

The random sub-path removal specifically includes: if there is an infeasible sub-path solution in the current solution, removing the infeasible sub-path; if the number of removed client nodes does not exceed the removal upper limit, selecting a random sub-path for removal. Remove until the number of removed client nodes reaches the removal limit;

The removal of the similar client nodes specifically includes: randomly selecting a client node i in the current solution (if there is an infeasible sub-path solution in the current solution, the infeasible sub-path solution is preferentially selected), and calculating the client node i Similarity with other customer nodes j, remove the customer nodes with the largest similarity in turn, until the number of removed customer nodes reaches the upper limit of removal. The distance d _ij between two customer points i and j is used as the basis for removing similar points. The smaller the value of d _ij is, the higher the similarity between customer points i and j is.

The removal of the similarity of the ring area specifically includes: considering the spatial relationship between the sub-paths in the current solution, and calculating the similarity based on three spatial relationships between the ring-shaped areas where each sub-path is located. The three spatial relationships are overlapping, intersecting, and separating. Figure 3 is a schematic diagram of the spatial relationship between the two sub-paths in the solution corresponding to the annular regions s ₁ and S ₂ . Each annular area is the smallest annular area that can cover the complete sub-path with the distribution center as the center, and the shaded area is the overlapping part of the corresponding annular areas of the two sub-paths.

Based on the three spatial relationships in Figure 3, the similarity is calculated;

The similarity between the two sub-paths is:

Among them, σ ₁ >σ ₂ >σ ₃ (its value is a fixed constant during the execution of the algorithm, which is a fixed value predetermined by humans. The selection of its value only needs to follow this inequality, which mainly represents a mutual relationship, and there is no fixed range standard , in this embodiment, σ ₁ =20, σ ₂ =5, σ ₃ =1), s _i and _{s j} represent the annular area corresponding to the two sub-paths ri and r _j in the current solution, then a single sub-path r The similarity coefficient of _i is:

中选出相似度系数最大的子路径r_m镠x，通过轮盘赌的方法选出非最大相似度系数的子路径r_d，之后移除r_m镠x和r_d中的所有客户节点，

表示所有子路径的相似度系数集合，l表示当前解中子路径总数。Using an elite strategy to remove subpaths,

Select the sub-path r _{m x x} with the largest similarity coefficient, select the sub-path r _d with the non-maximum similarity coefficient through the roulette method, and then remove all customer nodes in r _{m x} and r _d ,

Represents the set of similarity coefficients of all subpaths, and l represents the total number of subpaths in the current solution.

The cost-greedy insertion specifically includes: inserting each removed client node into the best position in the sub-path of the current solution. After inserting the customer node i to be inserted into the customer node j in the sub-path, its insertion cost c _ij =d _i,j +d _j,i+1 -d _i,i+1 , and the minimum insertion cost is selected for each insertion and subpaths that satisfy the load limit. If there are uninserted points after the above steps, only the insertion cost is considered for insertion until there are no remaining points.

The regret value insertion specifically includes: inserting each removed client node into the sub-path of the current solution, and selecting the sub-path with the largest regret value for each insertion; the regret value is: Δc=c _m2 -c _m1 , where c _m1 represents the smallest insertion cost and c _m2 represents the next smallest insertion cost. The significance is that the larger the regret value Δc, the greater the increase in the insertion cost caused by the selection of the sub-optimal position if the client node cannot be inserted into the optimal position.

A2. If the target value of the local solution is less than the target value of the current solution, the local solution is taken as the current solution; if the target value of the local solution is greater than the target value of the current solution, then according to the Metropolis criterion, combined with the temperature value of this iteration , the probability accepts the local solution and updates the temperature value; the temperature value in the first iteration process is the initial temperature value.

其中，P_Δ的值与[0,1]间的随机数r，若P_Δ≥r，则接受劣解作为下一次迭代的当前解，T表示当前温度，f(x)表示当前解x的目标值，f(x′)表示局部解x′的目标值。随机数r由程序随机函数生成。Probability Acceptance The acceptance probability of a local solution is

Among them, the value of P _Δ and the random number r between [0,1], if P _Δ ≥ r, the inferior solution is accepted as the current solution for the next iteration, T represents the current temperature, and f(x) represents the current solution x target value, f(x') represents the target value of the local solution x'. The random number r is generated by the program random function.

The update formula of temperature value is: T′=T*α

Among them, T′ represents the temperature after updating, T represents the current temperature before updating, α∈(0,1) is the cooling speed, and the value of α is a preset constant parameter, and the value is generally close to 1 (such as 0.95, 0.99 , 0.999, etc., here can be made 0.99), the closer to 1, the slower the temperature drop.

A3. According to the usage times and preset scores of each insertion operator and each removal operator, update the corresponding operator weight, and the update interval is the preset number of iterations. The preset number of iterations is generally artificially specified, and there is no hard requirement. In this embodiment, it is set to 20, which means that every iteration is 20 times, and the weight of each operator is adjusted based on the performance of each operator in these 20 iterations. The adjusted weights are used as the basis for the selection of the operator for the next 20 iterations.

The update formula of operator weight is:

表示算子q在τ次迭代内的使用次数，对应的预设评分为

在第y个迭代间隔内，第q个算子对应的算子权重为

在第y+1个迭代间隔内，第q个算子对应的算子权重为

μ∈[0,1]是控制算子影响的因子，该值越小，则第y个迭代间隔的权重对第y+1个迭代间隔的权重更新影响越大；若μ＝0，则权重完全取决于上一迭代间隔的权重；若μ＝1，则只考虑迭代间隔内评分，若0<μ<1，则同时考虑迭代间隔内的得分与之前的权重值，μ的值是人为预设的。where τ represents the interval τ iterations, y represents the y-th iteration interval, q represents the q-th operator,

represents the number of times the operator q is used in τ iterations, and the corresponding preset score is

In the yth iteration interval, the operator weight corresponding to the qth operator is

In the y+1th iteration interval, the operator weight corresponding to the qth operator is

μ∈[0,1] is the factor that controls the influence of the operator. The smaller the value, the greater the influence of the weight of the y-th iteration interval on the weight update of the y+1-th iteration interval; if μ=0, the weight It completely depends on the weight of the previous iteration interval; if μ=1, only the score in the iteration interval is considered, if 0<μ<1, the score in the iteration interval and the previous weight value are considered at the same time, and the value of μ is an artificial prediction. set.

The preset scoring rule of the operator is: if a new global optimal solution is found, the operator score is increased by θ ₁ ; if the new solution is only better than the current solution and worse than the global optimal solution, the operator score is increased by θ ₂ ; If the new solution is inferior to the current solution but is accepted by the Metropolis criterion as the initial solution for the next iteration, the operator score increases by θ ₃ ; if the new solution is inferior to the current solution and is not accepted by the Metropolis criterion, the operator score increases by θ ₄ . Generally, θ ₁ >θ ₂ >θ ₃ >θ ₄ , which is preset by humans and has no fixed standard. In this embodiment, θ ₁ =5, θ ₂ =3, θ ₃ =1.5, and θ ₄ =1.8.

A4. Stop the iteration until the target value of the local solution is not less than the current solution and the number of iterations reaches the preset number of times, and output the optimal coding sequence of the vehicle. Decoding the vehicle optimal code sequence is the vehicle optimal delivery route.

The above-mentioned preset number of times can be selected as 2000 generations, and the number of iterations of the unimproved solution can also be adjusted according to requirements. In addition, by selecting the recorded optimal solution for decoding, the vehicles required to complete all delivery tasks and the customers served by each vehicle can be obtained. As shown in Figure 4, Figure 4 shows the process of four vehicles to complete all delivery tasks and the customer nodes that each vehicle needs to serve.

This embodiment provides a method for planning a vehicle distribution path, including: acquiring basic information of a vehicle path; the basic information includes the geographical location of the distribution center node and each customer node, the number of customers, the demand for goods by customers, and the maximum load of the vehicle ; establish an objective function based on the distribution model; obtain an initial solution of the model through a scanning algorithm; take the objective function as the target value, and use an improved large-scale neighborhood search algorithm to select the optimal vehicle distribution path. Using the improved large-scale neighborhood search algorithm can filter out the optimal delivery route for vehicles, so as to formulate a more reasonable delivery route scheme, improve the delivery efficiency, and greatly reduce the total delivery mileage of vehicles. The improved large-scale neighborhood search algorithm in this embodiment takes into account the geographic location correlation of a single delivery in the delivery scheme, improves the local search capability of the algorithm, and has better applicability to the formulation of large-scale vehicle delivery schemes , which can provide a solution to the general vehicle distribution routing problem with capacity constraints, and provide a reference for managers to make decisions.

Embodiment 2:

As shown in FIG. 5 , this embodiment provides a vehicle distribution path planning system, the system includes:

The information acquisition module M1 is used to acquire the node location information, customer information and vehicle information of the distribution center; the customer information includes: the number of customers, the demand for goods of each customer, and the location of each customer node; the vehicle information includes: the number of vehicles and vehicle maximum load information;

The objective function determination module M2 is used to establish an objective function with the shortest total distribution vehicle mileage as the objective according to the distribution center node location information, customer information and vehicle maximum load information; the constraint condition of the objective function is that each customer consists of and be served by only one vehicle and the total customer demand served by each vehicle does not exceed the vehicle load;

The initial solution building module M3 is used to determine the initial coding sequence; the initial coding sequence is the coded representation of the initial distribution path, and the initial coding sequence includes each customer node on the initial distribution path and the order of delivery;

Initial solution target value calculation module M4, for calculating the target value of the initial coding sequence according to the target function;

The optimal distribution path iteration module M5 is used to calculate the initial temperature of the Metropolis criterion according to the target value of the initial coding sequence, and start the optimal distribution path iteration;

The optimal delivery route iteration specifically includes:

According to the weight of each removal operator and each insertion operator, use the roulette method to select a removal operator and an insertion operator, and according to the selected removal operator and insertion operator, respectively, for the current solution Perform an operator removal operation and an operator insertion operation to obtain a local solution; the removal operators include: random customer node removal, random subpath removal, similar customer node removal and ring similarity removal; all The insertion operator includes: cost greedy insertion and regret value insertion; the current solution in the first iteration process is the initial coding sequence;

If the target value of the local solution is smaller than the target value of the current solution, the local solution is taken as the current solution; if the target value of the local solution is greater than the target value of the current solution, then according to the Metropolis criterion, combined with the temperature value of this iteration, the probability Accept the local solution and update the temperature value; the temperature value in the first iteration process is the initial temperature value;

According to the usage times and preset scores of each insertion operator and each removal operator, the corresponding operator weights are updated, and the update interval is the preset number of iterations;

Until the target value of the local solution is not less than the current solution iteration times and reaches the preset number of times, the iteration is stopped, and the vehicle optimal coding sequence is output.

The various embodiments in this specification are described in a progressive manner, and each embodiment focuses on the differences from other embodiments, and the same and similar parts between the various embodiments can be referred to each other. For the system disclosed in the embodiment, since it corresponds to the method disclosed in the embodiment, the description is relatively simple, and the relevant part can be referred to the description of the method.

In this paper, specific examples are used to illustrate the principles and implementations of the present invention. The descriptions of the above embodiments are only used to help understand the methods and core ideas of the present invention; meanwhile, for those skilled in the art, according to the present invention There will be changes in the specific implementation and application scope. In conclusion, the contents of this specification should not be construed as limiting the present invention.

Claims (10)

1. A vehicle distribution path planning method, characterized in that the method comprises: Acquiring the location information, customer information and vehicle information of the distribution center node; the customer information includes: the number of customers, the demand for goods of each customer, and the location of each customer node; the vehicle information includes: the number of vehicles and the maximum load information of the vehicle; According to the distribution center node location information, customer information and vehicle maximum load information, establish an objective function with the shortest total distribution vehicle mileage as the goal; the constraint condition of the objective function is that each customer is served by one and only one vehicle and each The total customer demand served by the vehicle does not exceed the vehicle load; determining an initial coding sequence; the initial coding sequence is a coded representation of the initial delivery route, and the initial coding sequence includes each client node on the initial delivery route and the delivery sequence of each client node; Calculate the target value of the initial coding sequence according to the target function; Calculate the initial temperature of the Metropolis criterion according to the target value of the initial coding sequence, and start the optimal delivery route iteration; The optimal delivery route iteration specifically includes: According to the weight of each removal operator and each insertion operator, use the roulette method to select a removal operator and an insertion operator, and according to the selected removal operator and insertion operator, respectively, for the current solution Perform an operator removal operation and an operator insertion operation to obtain a local solution; the removal operators include: random customer node removal, random subpath removal, similar customer node removal and ring similarity removal; all The insertion operator includes: cost greedy insertion and regret value insertion; the current solution in the first iteration process is the initial coding sequence; If the target value of the local solution is smaller than the target value of the current solution, the local solution is taken as the current solution; if the target value of the local solution is greater than the target value of the current solution, then according to the Metropolis criterion, combined with the temperature value of this iteration, the probability Accept the local solution and update the temperature value; the temperature value in the first iteration process is the initial temperature value; According to the usage times and preset scores of each insertion operator and each removal operator, the corresponding operator weights are updated, and the update interval is the preset number of iterations; Until the target value of the local solution is not less than the current solution iteration times and reaches the preset number of times, the iteration is stopped, and the vehicle optimal coding sequence is output. 2. The vehicle distribution route planning method according to claim 1, characterized in that, the objective function that takes the shortest total distribution vehicle mileage as the goal is:

Among them, F represents the objective function; m is the maximum number of vehicles that can be used for delivery; n is the number of customers; _{d ij} represents the distance between any two nodes i and j; i and j, 0 means that vehicle k does not pass through nodes i and j, and 1 means that vehicle k passes through nodes i and j. 3. The vehicle distribution path planning method according to claim 1, wherein the determining the initial coding sequence specifically comprises: Taking the position of the distribution center node as the pole, and taking the ray of the distribution center node and the customer node closest to the distribution center node as the polar axis, a polar coordinate system is established, and the position coordinates of the distribution center node and the customer node are converted. is polar coordinates; Sort client nodes in order of polar angle size; On the premise that the constraints are met, each customer node is inserted into the solution in turn according to the objective function; the solution includes several distribution paths; during the insertion process, if a customer node that does not meet the constraints is encountered, and all If the client node that does not meet the constraints is not the last client node, the client node that does not meet the constraints is skipped; Search whether there is an uninserted customer node j, and if so, calculate the insertion cost c _ij between each two customer nodes inserted by the customer node j into each delivery route, and insert the two corresponding to the lowest cost. between client nodes; where c _ij =d _i,j +d _j,i+1 -d _i,i+1 ; in the formula, d _i,j represents the distance between client node i and client node j, d _j,i+1 represents the distance between client node j and client node i+1, and d _i,i+1 represents the distance between client node i and client node i+1; Traverse all uninserted client nodes until all client nodes are inserted into the solution to obtain the initial encoding sequence. 4. The vehicle distribution path planning method according to claim 1, wherein the calculating the initial temperature of the Metropolis criterion according to the target value in the initial coding sequence, specifically comprises: Calculate the initial temperature T ₀ of the Metropolis criterion according to the target value in the initial coding sequence;

The initial temperature T ₀ indicates that under the condition that the acceptance probability of the inferior solution is η, the inferior solution that is worse than the current target value f(x) by β% is to be accepted, and the corresponding temperature is T ₀ . 5. The vehicle distribution path planning method according to claim 1, wherein, The random customer node removal specifically includes: if there is an infeasible sub-path solution in the current solution, randomly removing the customer nodes in the infeasible sub-path until the sum of the customer demand does not exceed the maximum vehicle load; If the number of customer nodes does not exceed the removal upper limit, select random customer nodes for removal until the number of removed customer nodes reaches the removal upper limit; the solution of the infeasible sub-path is that the sum of the customer demands of the sub-paths is greater than the maximum vehicle load; the removal cap is generated by a random function; The random sub-path removal specifically includes: if there is an infeasible sub-path solution in the current solution, removing the infeasible sub-path; if the number of removed client nodes does not exceed the removal upper limit, selecting a random sub-path for removal. Remove until the number of removed client nodes reaches the removal limit; The removal of similar client nodes specifically includes: randomly selecting a client node i in the current solution, calculating the similarity between the client node i and other client nodes j, and sequentially removing the client node with the largest similarity until the client node is removed. The number of nodes reaches the removal limit; The removal of the similarity of the ring area specifically includes: calculating the similarity based on three spatial relationships between the ring areas where each sub-path in the current solution is located; The similarity between the two sub-paths is:

Among them, σ ₁ >σ ₂ >σ ₃ , s _i and _{s j} represent the annular area corresponding to the two sub-paths ri and r _j in the current solution, then the similarity coefficient of a single sub-path _ri is:

exist

The sub-path r _max with the largest similarity coefficient is selected from among them, and the sub-path r _d with the non-maximum similarity coefficient is selected by the roulette method, and then all customer nodes in r _{m, x} and r _d are removed,

Represents the set of similarity coefficients of all subpaths, and l represents the total number of subpaths in the current solution; The cost-greedy insertion specifically includes: inserting each removed client node into the sub-path of the current solution, and selecting the sub-path with the smallest insertion cost and satisfying the load restriction for each insertion; The regret value insertion specifically includes: inserting each removed client node into the sub-path of the current solution, and selecting the sub-path with the largest regret value for each insertion; the regret value is: Δc=c _m2 -c _m1 , where c _m1 represents the smallest insertion cost and c _m2 represents the next smallest insertion cost. 6. The vehicle distribution path planning method according to claim 1, wherein the update formula of the operator weight is:

where τ represents the interval τ iterations, y represents the y-th iteration interval, q represents the q-th operator,

represents the number of times the operator q is used in τ iterations, and the corresponding preset score is

In the yth iteration interval, the operator weight corresponding to the qth operator is

In the y+1th iteration interval, the operator weight corresponding to the qth operator is

μ∈[0,1] is the factor that controls the influence of the operator. The smaller the value, the greater the influence of the weight of the y-th iteration interval on the weight update of the y+1-th iteration interval; if μ=0, the weight It depends entirely on the weight of the previous iteration interval; if μ=1, only the score in the iteration interval is considered, and if 0<μ<1, the score in the iteration interval and the previous weight value are considered at the same time. 7. The vehicle distribution path planning method according to claim 1, wherein before said determining the initial coding sequence, further comprising: Chromosome coding is performed on the distribution center node and each of the customer nodes. 8 . The vehicle distribution route planning method according to claim 1 , wherein the update formula of the temperature value is: T′=T*α Among them, T′ represents the temperature after updating, T represents the current temperature before updating, α∈(0,1) is the cooling speed, and the value of α is a preset constant parameter. 9 . The vehicle distribution route planning method according to claim 1 , wherein, according to the Metropolis criterion, combined with the temperature value of this iteration, the probability of accepting a local solution includes: According to the Metropolis criterion, combined with the temperature value of this iteration, the probability of accepting the local solution is

Among them, the value of P _Δ and the random number r between [0,1], if P _Δ ≥ r, the inferior solution is accepted as the current solution for the next iteration, T represents the current temperature, and f(x) represents the current solution x target value, f(x') represents the target value of the local solution x'. 10. A vehicle distribution path planning system, wherein the system comprises: The information acquisition module is used to acquire the node location information, customer information and vehicle information of the distribution center; the customer information includes: the number of customers, the demand for goods of each customer, and the location of each customer node; the vehicle information includes: the number of vehicles and Vehicle maximum load information; The objective function determination module is used to establish an objective function with the shortest total distribution vehicle mileage as the objective according to the distribution center node location information, customer information and vehicle maximum load information; the constraint condition of the objective function is that each customer consists of and only be served by one vehicle and the total customer demand served by each vehicle does not exceed the vehicle load; The initial solution building module is used to determine the initial coding sequence; the initial coding sequence is the coded representation of the initial distribution path, and the initial coding sequence includes each customer node on the initial distribution path and the distribution of each customer node order; an initial solution target value calculation module, for calculating the target value of the initial coding sequence according to the target function; an optimal delivery route iteration module, configured to calculate the initial temperature of the Metropolis criterion according to the target value of the initial coding sequence, and start the optimal delivery route iteration; The optimal delivery route iteration specifically includes: According to the weight of each removal operator and each insertion operator, use the roulette method to select a removal operator and an insertion operator, and according to the selected removal operator and insertion operator, respectively, for the current solution Perform operator removal operation and operator insertion operation to obtain a local solution; the removal operator includes: random customer node removal, random subpath removal, similar customer node removal and loop similarity removal; all The insertion operator includes: cost greedy insertion and regret value insertion; the current solution in the first iteration process is the initial coding sequence; If the target value of the local solution is smaller than the target value of the current solution, the local solution is taken as the current solution; if the target value of the local solution is greater than the target value of the current solution, then according to the Metropolis criterion, combined with the temperature value of this iteration, the probability Accept the local solution and update the temperature value; the temperature value in the first iteration process is the initial temperature value; According to the usage times and preset scores of each insertion operator and each removal operator, update the corresponding operator weight, and the update interval is the preset number of iterations; Until the target value of the local solution is not less than the current solution iteration times and reaches the preset number of times, the iteration is stopped, and the vehicle optimal coding sequence is output.

Images