CN103530699A - Multi-time-window vehicle route selection method on basis of improved universal gravitation algorithm - Google Patents

Abstract

Translated fromChinese

本发明公开了一种基于改进万有引力算法的多时间窗车辆路径选择方法，主要解决现有技术不能同时对离散型和聚集型多时间窗车辆路径进行选择的问题。其实现步骤为：（1）对输入的客户进行聚类处理，确定聚类中心，并对车辆编号；（2）以车场为中心，根据聚类中心，将客户分布区域分为车辆数个片区，确定每辆车服务的客户集合；（3）记第k辆车服务的客户集为一个种群，在t＝0时刻，初始化该种群，并用万有引力算法对该种群进行寻优；（4）对该种群中各个体执行交叉操作，得到一个临时种群，再对该临时种群的各个体执行边界约束检查，并通过下、上层终止判断，得到多时间窗车辆路径问题的最优路径。本发明能同时对不同类型的多时间窗车辆路径进行选择。

The invention discloses a multi-time-window vehicle route selection method based on an improved gravitational algorithm, which mainly solves the problem that the prior art cannot simultaneously select discrete and aggregated multi-time-window vehicle routes. The implementation steps are: (1) cluster the input customers, determine the cluster center, and number the vehicles; (2) take the parking lot as the center, divide the customer distribution area into several vehicle areas according to the cluster center , determine the set of customers served by each car; (3) record the set of customers served by the kth car as a population, initialize the population at time t=0, and use the gravitational algorithm to optimize the population; (4) Each individual in the population performs a crossover operation to obtain a temporary population, and then performs a boundary constraint check on each individual of the temporary population, and judges the termination of the lower and upper layers to obtain the optimal path of the multi-time window vehicle routing problem. The invention can simultaneously select different types of multi-time window vehicle routes.

Description

Translated fromChinese

基于改进万有引力算法的多时间窗车辆路径选择方法Multi-time-window vehicle routing method based on improved gravitational algorithm

技术领域technical field

本发明属于交通运输技术领域，特别涉及多时间窗的车辆路径选择方法，可用于对离散型和聚集型的多时间窗车辆路径进行调度。The invention belongs to the technical field of transportation, in particular to a multi-time-window vehicle route selection method, which can be used to schedule discrete and aggregate multi-time-window vehicle routes.

背景技术Background technique

多时间窗车辆路径问题是指车辆从配送中心出发服务客户，需要在用户提供若干个时间窗内，选择在唯一的一个时间窗内抵达服务，要求每个客户只能被一辆车服务且仅服务一次，路径选择目标是使得在满足用户时间要求及车辆载荷有限等条件下成本最小。基于多时间窗的车辆路径问题广泛存在于当今现实生活的物流运输之中，而国内对于此问题的研究还较少且一般仅只能解决某种特殊类型的多时间窗车辆调度问题。The multi-time-window vehicle routing problem means that the vehicle departs from the distribution center to serve customers. It needs to choose to arrive at the service in a unique time window within several time windows provided by the user. It is required that each customer can only be served by one vehicle and only Service once, the goal of route selection is to minimize the cost under the conditions of meeting the user's time requirements and limited vehicle load. The vehicle routing problem based on multiple time windows widely exists in today's real-life logistics transportation, but domestic research on this problem is still less and generally only solves a special type of multi-time window vehicle scheduling problem.

多时间窗车辆路径问题是一个NP-Hard问题，这意味着在问题规模增大到一定数量时将很难或者根本无法求得问题的全局最优解。采用精确算法虽然可以对小规模的多时间窗车辆调度问题得到最优解，但却不适用于求解现实中的大规模的多时间窗车辆调度问题。有些学者运用序列插入启发式算法求解多时间窗车辆路径问题，虽然取得不错结果，但此种算法只能对聚集型或离散型多时间窗车辆路径选择起作用，而不能同时适用于两种类型的多时间窗车辆路径选择。The multi-time window vehicle routing problem is an NP-hard problem, which means that when the problem scale increases to a certain number, it will be difficult or impossible to obtain the global optimal solution of the problem. Although the exact algorithm can obtain the optimal solution to the small-scale multi-time window vehicle scheduling problem, it is not suitable for solving the large-scale multi-time window vehicle scheduling problem in reality. Some scholars use the sequence insertion heuristic algorithm to solve the multi-time window vehicle routing problem. Although good results have been achieved, this algorithm can only work for aggregated or discrete multi-time window vehicle routing, and cannot be applied to both types at the same time. Multi-time window vehicle routing.

发明内容Contents of the invention

本发明的目的在于针对上述已有技术的不足，提出一种基于万有引力算法的多时间窗车辆路径选择方法，以同时对两种类型的多时间窗车辆路径进行选择。The object of the present invention is to address the above-mentioned deficiencies in the prior art, and propose a multi-time-window vehicle route selection method based on a gravitational algorithm, so as to simultaneously select two types of multi-time-window vehicle routes.

实现本发明目的技术思路是对现有的万有引力算法GSA进行改进，以有效解决离散型和聚集型多时间窗的车辆路径问题，提高搜索能力，其技术方案包括如下步骤：The technical idea of realizing the object of the present invention is to improve the existing gravitational algorithm GSA, to effectively solve the vehicle path problem of discrete and aggregated multi-time windows, and improve the search ability, and its technical solution includes the following steps:

（1）输入客户分布位置情况,车场中的车辆数M，下层迭代次数D和上层调整次数U，对客户进行快速聚类处理，即随机选择M个点作为聚类中心，设聚类调整次数I＝1；(1) Input the location of customers, the number of vehicles M in the yard, the number of iterations of the lower layer D and the number of adjustments U of the upper layer, and quickly cluster the customers, that is, randomly select M points as the cluster centers, and set the number of cluster adjustments I=1;

（2）聚类分组：对每辆车进行编号，令第一辆车的序号k＝1；以车场为中心，根据聚类中心，将客户分布区域分为车辆数个片区，确定每辆车服务的客户集合；(2) Clustering and grouping: Number each vehicle, so that the serial number of the first vehicle is k=1; take the parking lot as the center, divide the customer distribution area into several vehicle areas according to the clustering center, and determine each vehicle the collection of clients served;

（3）初始化种群：记第k辆车服务的客户集为一个种群，在t＝0时刻，初始化该种群中每个个体的空间位置；(3) Initialize the population: record the customer set served by the kth car as a population, and initialize the spatial position of each individual in the population at time t=0;

（4）利用万有引力算法寻优：在t时刻，选取种群中适应度值最好的m个个体，经万有引力作用使其施力于种群中的其他个体，然后更新一次个体速度和个体位置，得到T＝t+1时刻种群初始速度和位置，令t＝t+1；(4) Use the gravitational algorithm to optimize: at time t, select the m individuals with the best fitness value in the population, and make them exert force on other individuals in the population through the action of gravity, and then update the individual speed and individual position once, and get The initial velocity and position of the population at T=t+1, let t=t+1;

（5）局部搜索：(5) Local search:

（5a）对更新后的种群中每个个体进行一次交叉操作，得到一个临时种群，再对临时种群中每个个体执行边界约束检查，若该种群中存在目标值优于初始个体目标值的个体，则用其替换初始个体，否则种群中每个个体的位置不变；(5a) Perform a crossover operation on each individual in the updated population to obtain a temporary population, and then perform a boundary constraint check on each individual in the temporary population, if there is an individual in the population whose target value is better than the initial individual target value , then use it to replace the initial individual, otherwise the position of each individual in the population remains unchanged;

（5b）比较T与下层迭代次数D，若T≤D,返回步骤（4）；否则，记录第k辆车的最优路径、路径长及迟到时间；(5b) Compare T with the number of iterations D of the lower layer, if T≤D, return to step (4); otherwise, record the optimal path, path length and late time of the kth car;

（6）下层终止判断：将车辆序号k+1与车辆数M进行比较，若(k+1)≤M，则返回步骤（3）；否则，计算当前分组下寻得的最优个体，并将该最优个体与种群中已知最优个体相比较，若当前分组下寻得的最优个体的目标函数值小于种群中已知最优个体的目标函数值，则用当前分组下寻得的最优个体替换种群中已知最优个体；否则，保持种群中已知最优个体不变；(6) Termination judgment of the lower layer: compare the vehicle serial number k+1 with the number of vehicles M, if (k+1)≤M, return to step (3); otherwise, calculate the optimal individual found under the current group, and Compare the optimal individual with the known optimal individual in the population, if the objective function value of the optimal individual found under the current group is smaller than the objective function value of the known optimal individual in the population, then use the optimal individual found under the current group Replace the known best individual in the population with the optimal individual; otherwise, keep the known optimal individual in the population unchanged;

（7）上层终止判断：将聚类调整次数I+1与上层调整次数U进行比较，若(I+1)≤U，则对输入的客户重新进行聚类，确定聚类中心，并返回步骤（2）；否则，停止聚类，并返回记录的已知最优路径、路径长及迟到时间，该路径即为所选多时间窗车辆路径的最优路径。(7) Upper-level termination judgment: compare the number of clustering adjustments I+1 with the number of upper-level adjustments U, if (I+1)≤U, re-cluster the input customers, determine the clustering center, and return to the step (2); otherwise, stop clustering and return the recorded known optimal route, route length and late time, which is the optimal route of the selected multi-time window vehicle route.

本发明与现有技术相比具有以下优点：Compared with the prior art, the present invention has the following advantages:

1、本发明通过聚类和改进的万有引力算法，能实现对离散型和聚集型的多时间窗车辆路径选择，而目前国内对于多时间窗车辆路径选择问题的研究较少，且一般研究只能解决特定的一种类型的多时间窗车辆路径问题。1. Through clustering and improved gravitational algorithm, the present invention can realize the selection of discrete and aggregated multi-time window vehicle routes. However, there are few domestic studies on the problem of multi-time window vehicle route selection, and the general research can only Solve a specific type of vehicle routing problem with multiple time windows.

2、本发明与现有的万有引力算法相比较，表现出更强的搜索能力，而且有效性和实用性突出。2. Compared with the existing gravitational algorithm, the present invention shows stronger search ability, and has outstanding effectiveness and practicability.

附图说明Description of drawings

图1是本发明的实现步骤流程图。Fig. 1 is a flowchart of the implementation steps of the present invention.

具体实施方法：Specific implementation method:

参照图1，本发明的实现步骤如下：With reference to Fig. 1, the realization steps of the present invention are as follows:

步骤1、输入客户分布位置情况，车场中的车辆总数M和客户总数C，设下层迭代次数为D≥5，上层调整次数为U≥5，聚类调整次数为I＝1，对客户进行快速聚类处理，即随机选择M个点作为聚类中心(X_i,Y_i),i＝1,2,…,M。Step 1. Input the location of customers, the total number of vehicles M and the total number of customers C in the yard, set the number of iterations in the lower layer to be D≥5, the number of adjustments in the upper layer to be U≥5, and the number of clustering adjustments to be I=1. Clustering processing, that is, randomly select M points as cluster centers (X_i , Y_i ), i=1, 2, . . . , M.

步骤2、聚类分组Step 2, clustering and grouping

对每辆车进行编号，令第一辆车的序号k＝1；以车场为中心，根据聚类中心(X_i,Y_i),i＝1,2,…,M，用下式将客户分布区域分为M个片区，确定每辆车服务的客户集合：Number each car, make the serial number of the first car k=1; take the parking lot as the center, according to the cluster center (X_i ,Y_i ), i=1,2,...,M, use the following formula to classify customers The distribution area is divided into M areas, and the set of customers served by each vehicle is determined:

$\mathrm{<span><span>min</span>min</span>}\underset{\mathrm{<span><span>r</span>r</span>}<span><span>=</span>=</span>\mathrm{<span><span>1</span>1</span>}}{\overset{\mathrm{<span><span>M</span>m</span>}}{\mathrm{<span><span>\&Sigma;</span>\&Sigma;</span>}}}\underset{\mathrm{<span><span>i</span>i</span>}<span><span>\&Element;</span>\&Element;</span>{\mathrm{<span><span>S</span>S</span>}}_{\mathrm{<span><span>r</span>r</span>}}}{\mathrm{<span><span>\&Sigma;</span>\&Sigma;</span>}}\sqrt{{<span><span>(</span>(</span>{\mathrm{<span><span>x</span>x</span>}}_{\mathrm{<span><span>i</span>i</span>}}<span><span>-</span>-</span>{\mathrm{<span><span>X</span>x</span>}}_{\mathrm{<span><span>r</span>r</span>}}<span><span>)</span>)</span>}^{\mathrm{<span><span>2</span>2</span>}}<span><span>+</span>+</span>{<span><span>(</span>(</span>{\mathrm{<span><span>y</span>the\; y</span>}}_{\mathrm{<span><span>i</span>i</span>}}<span><span>-</span>-</span>{\mathrm{<span><span>Y</span>Y</span>}}_{\mathrm{<span><span>r</span>r</span>}}<span><span>)</span>)</span>}^{\mathrm{<span><span>2</span>2</span>}}}<span><span>,</span>,</span>$

s.t.S_i∩S_j＝φ,

j∈{1,2,…,M}，且i≠j，stS_i ∩S_j = φ,

j∈{1,2,…,M}, and i≠j,

其中(x_i,y_i)表示第i个客户的位置，S_r表示第r辆车服务的客户集合，C表示需要服务的客户总数，s.t.表示受约束，∩表示集合的交，

表示任意，∪表示集合的并。Where (_xi , y_i ) represents the position of the i-th customer, S_r represents the set of customers served by the r-th car, C represents the total number of customers who need to be served, st represents constrained, ∩ represents the intersection of the sets,

means any, and ∪ means the union of sets.

步骤3、初始化种群Step 3. Initialize the population

记第k辆车服务的客户集为一个种群，在t＝0时刻，用以下公式对该种群中每个个体进行初始化：Record the customer set of the kth car service as a population, and at time t=0, use the following formula to initialize each individual in the population:

${\mathrm{<span><span>p</span>p</span>}}_{\mathrm{<span><span>j</span>j</span>}}^{\mathrm{<span><span>i</span>i</span>}}<span><span>=</span>=</span><span><span>(</span>(</span>\mathrm{<span><span>b</span>b</span>}<span><span>-</span>-</span>\mathrm{<span><span>a</span>a</span>}<span><span>)</span>)</span>\mathrm{<span><span>R</span>R</span>}<span><span>+</span>+</span>\mathrm{<span><span>a</span>a</span>}<span><span>,</span>,</span>{\mathrm{<span><span>p</span>p</span>}}_{\mathrm{<span><span>j</span>j</span>}}^{\mathrm{<span><span>i</span>i</span>}}<span><span>\&Element;</span>\&Element;</span><span><span>[</span>[</span>\mathrm{<span><span>a</span>a</span>}<span><span>,</span>,</span>\mathrm{<span><span>b</span>b</span>}<span><span>]</span>]</span><span><span>,</span>,</span>\mathrm{<span><span>i</span>i</span>}<span><span>=</span>=</span>\mathrm{<span><span>1,2</span>1,2</span>}<span><span>,</span>,</span><span><span>.</span>.</span><span><span>.</span>.</span><span><span>.</span>.</span><span><span>,</span>,</span>\mathrm{<span><span>L</span>L</span>}<span><span>,</span>,</span>\mathrm{<span><span>j</span>j</span>}<span><span>=</span>=</span>\mathrm{<span><span>1,2</span>1,2</span>}<span><span>,</span>,</span><span><span>.</span>.</span><span><span>.</span>.</span><span><span>.</span>.</span><span><span>,</span>,</span>{\mathrm{<span><span>N</span>N</span>}}_{\mathrm{<span><span>k</span>k</span>}}$

其中，表示第i个个体在第j维空间的位置坐标，R表示[0,1]内的随机数，L表示种群规模，N_k表示第k辆车服务客户数目，[a,b]表示空间范围。in, Represents the position coordinates of the i-th individual in the j-dimensional space, R represents a random number in [0,1], L represents the population size, N_k represents the number of customers served by the k-th car, [a, b] represents the spatial range .

步骤4、利用万有引力算法寻优Step 4. Use the gravitational algorithm to optimize

（4a）在t时刻，定义第i个个体和第j个个体在第d维空间的万有引力为：(4a) At time t, define the gravitational force of the i-th individual and the j-th individual in the d-dimensional space as:

${\mathrm{<span><span>F</span>f</span>}}_{\mathrm{<span><span>ij</span>ij</span>}}^{\mathrm{<span><span>d</span>d</span>}}<span><span>(</span>(</span>\mathrm{<span><span>t</span>t</span>}<span><span>)</span>)</span><span><span>=</span>=</span>\mathrm{<span><span>G</span>G</span>}<span><span>(</span>(</span>\mathrm{<span><span>t</span>t</span>}<span><span>)</span>)</span>\frac{{\mathrm{<span><span>M</span>m</span>}}_{\mathrm{<span><span>pi</span>p</span>}}<span><span>(</span>(</span>\mathrm{<span><span>t</span>t</span>}<span><span>)</span>)</span><span><span>\&times;</span>\&times;</span>{\mathrm{<span><span>M</span>m</span>}}_{\mathrm{<span><span>oj</span>oj</span>}}<span><span>(</span>(</span>\mathrm{<span><span>t</span>t</span>}<span><span>)</span>)</span>}{{\mathrm{<span><span>R</span>R</span>}}_{\mathrm{<span><span>ij</span>ij</span>}}<span><span>(</span>(</span>\mathrm{<span><span>t</span>t</span>}<span><span>)</span>)</span><span><span>+</span>+</span>\mathrm{<span><span>\&epsiv;</span>\&epsiv;</span>}}<span><span>(</span>(</span>{\mathrm{<span><span>p</span>p</span>}}_{\mathrm{<span><span>j</span>j</span>}}^{\mathrm{<span><span>d</span>d</span>}}<span><span>(</span>(</span>\mathrm{<span><span>t</span>t</span>}<span><span>)</span>)</span><span><span>-</span>-</span>{\mathrm{<span><span>p</span>p</span>}}_{\mathrm{<span><span>i</span>i</span>}}^{\mathrm{<span><span>d</span>d</span>}}<span><span>(</span>(</span>\mathrm{<span><span>t</span>t</span>}<span><span>)</span>)</span><span><span>)</span>)</span><span><span>,</span>,</span>$

其中，G(t)是t时刻的万有引力系数，R_ij(t)表示t时刻第i个个体与第j个个体之间的欧几里德距离，ε是个充分小的常数，M_pi(t)表示t时刻受力第i个个体的质量，M_oj(t)表示t时刻施力第j个个体的质量，

(t)表示第i个个体在第d维空间中的位置；Among them, G(t) is the gravitational coefficient at time t, R_ij (t) represents the Euclidean distance between the i-th individual and the j-th individual at time t, ε is a sufficiently small constant, M_pi (t ) represents the mass of the i-th individual under force at time t, M_oj (t) represents the mass of the j-th individual under force at time t,

(t) represents the position of the i-th individual in the d-dimensional space;

（4b）根据上述万有引力公式，分别计算t时刻，第i个个体在第d维空间的合力F_i^d和第i个个体在第d维空间产生的加速度

(t)：(4b) According to the above formula of universal gravitation, respectively calculate the resultant force F_i^d of the i-th individual in the d-dimensional space and the acceleration produced by the i-th individual in the d-dimensional space at time t

(t):

${{\mathrm{<span><span>F</span>f</span>}}_{\mathrm{<span><span>i</span>i</span>}}}^{\mathrm{<span><span>d</span>d</span>}}<span><span>=</span>=</span>\underset{\mathrm{<span><span>j</span>j</span>}<span><span>=</span>=</span>\mathrm{<span><span>1</span>1</span>}<span><span>,</span>,</span>\mathrm{<span><span>j</span>j</span>}<span><span>\&NotEqual;</span>\&NotEqual;</span>\mathrm{<span><span>i</span>i</span>}}{\overset{\mathrm{<span><span>N</span>N</span>}}{\mathrm{<span><span>\&Sigma;</span>\&Sigma;</span>}}}\mathrm{<span><span>R</span>R</span>}<span><span>\&times;</span>\&times;</span>{\mathrm{<span><span>F</span>f</span>}}_{\mathrm{<span><span>ij</span>ij</span>}}^{\mathrm{<span><span>d</span>d</span>}}<span><span>(</span>(</span>\mathrm{<span><span>t</span>t</span>}<span><span>)</span>)</span>$

${\mathrm{<span><span>a</span>a</span>}}_{\mathrm{<span><span>i</span>i</span>}}^{\mathrm{<span><span>d</span>d</span>}}<span><span>(</span>(</span>\mathrm{<span><span>t</span>t</span>}<span><span>)</span>)</span><span><span>=</span>=</span>\frac{{\mathrm{<span><span>F</span>f</span>}}_{\mathrm{<span><span>i</span>i</span>}}^{\mathrm{<span><span>d</span>d</span>}}<span><span>(</span>(</span>\mathrm{<span><span>t</span>t</span>}<span><span>)</span>)</span>}{{\mathrm{<span><span>M</span>m</span>}}_{\mathrm{<span><span>ii</span>i</span>}}<span><span>(</span>(</span>\mathrm{<span><span>t</span>t</span>}<span><span>)</span>)</span>}<span><span>,</span>,</span>$

其中，R是[0,1]内的随机数，N表示在第d维空间中对第i个个体施加万有引力的个体数，(t)为t时刻第i个个体和第j个个体之间的万有引力，M_ii(t)表示t时刻第i个个体的质量；Among them, R is a random number in [0,1], and N represents the number of individuals exerting universal gravitation on the i-th individual in the d-th dimensional space, (t) is the gravitational force between the i-th individual and the j-th individual at time t, and M_ii (t) represents the mass of the i-th individual at time t;

（4c）根据上述加速度公式，更新个体速度和个体位置：(4c) According to the above acceleration formula, update the individual velocity and individual position:

在t时刻，选取种群中适应度值最好的m个个体，经万有引力作用使得其施力于种群中的其他个体，然后更新一次个体速度

(t)和个体位置

(t)，得到T＝t+1时刻种群中每个个体的速度

和位置

At time t, select the m individuals with the best fitness value in the population, make them exert force on other individuals in the population through the action of gravitation, and then update the individual speed

(t) and individual position

(t), get the speed of each individual in the population at time T=t+1

and location

${\mathrm{<span><span>v</span>v</span>}}_{\mathrm{<span><span>i</span>i</span>}}^{\mathrm{<span><span>d</span>d</span>}}<span><span>(</span>(</span>\mathrm{<span><span>t</span>t</span>}<span><span>+</span>+</span>\mathrm{<span><span>1</span>1</span>}<span><span>)</span>)</span><span><span>=</span>=</span>\mathrm{<span><span>R</span>R</span>}<span><span>\&times;</span>\&times;</span>{\mathrm{<span><span>v</span>v</span>}}_{\mathrm{<span><span>i</span>i</span>}}^{\mathrm{<span><span>d</span>d</span>}}<span><span>(</span>(</span>\mathrm{<span><span>t</span>t</span>}<span><span>)</span>)</span><span><span>+</span>+</span>{\mathrm{<span><span>a</span>a</span>}}_{\mathrm{<span><span>i</span>i</span>}}^{\mathrm{<span><span>d</span>d</span>}}<span><span>(</span>(</span>\mathrm{<span><span>t</span>t</span>}<span><span>)</span>)</span>$

${\mathrm{<span><span>p</span>p</span>}}_{\mathrm{<span><span>i</span>i</span>}}^{\mathrm{<span><span>d</span>d</span>}}<span><span>(</span>(</span>\mathrm{<span><span>t</span>t</span>}<span><span>+</span>+</span>\mathrm{<span><span>1</span>1</span>}<span><span>)</span>)</span><span><span>=</span>=</span>{\mathrm{<span><span>p</span>p</span>}}_{\mathrm{<span><span>i</span>i</span>}}^{\mathrm{<span><span>d</span>d</span>}}<span><span>(</span>(</span>\mathrm{<span><span>t</span>t</span>}<span><span>)</span>)</span><span><span>+</span>+</span>{\mathrm{<span><span>v</span>v</span>}}_{\mathrm{<span><span>i</span>i</span>}}^{\mathrm{<span><span>d</span>d</span>}}<span><span>(</span>(</span>\mathrm{<span><span>t</span>t</span>}<span><span>+</span>+</span>\mathrm{<span><span>1</span>1</span>}<span><span>)</span>)</span><span><span>,</span>,</span>$

其中

为t时刻第i个个体在第d维空间中的速度，

为t时刻第i个个体在第d维空间中的加速度，

为t时刻第i个个体在第d维空间中的位置，R是[0,1]内的随机数。in

is the velocity of the i-th individual in the d-dimensional space at time t,

is the acceleration of the i-th individual in the d-dimensional space at time t,

is the position of the i-th individual in the d-dimensional space at time t, and R is a random number in [0,1].

步骤5、局部搜索Step 5, local search

（5a）令t＝t+1，对更新后的种群中每个个体进行一次交叉操作，得到一个临时种群，再对临时种群中每个个体执行边界约束检查，若该种群中存在目标值优于初始个体目标值的个体，则用其替换初始个体，否则种群中每个个体的位置不变；(5a) Let t=t+1, perform a crossover operation on each individual in the updated population to obtain a temporary population, and then perform a boundary constraint check on each individual in the temporary population, if there is a target value superiority in the population If the individual is lower than the target value of the initial individual, it is used to replace the initial individual, otherwise the position of each individual in the population remains unchanged;

（5b）比较T与下层迭代次数D，若T≤D,返回步骤4；否则，记录第k辆车的最优路径、路径长及迟到时间。(5b) Compare T with the number of iterations D of the lower layer, if T≤D, return to step 4; otherwise, record the optimal path, path length and late time of the kth car.

步骤6、下层终止判断Step 6. Lower layer termination judgment

将车辆序号k+1与车辆数M进行比较，若(k+1)≤M，则令k＝k+1并返回步骤3；否则，计算当前分组下寻得的最优个体，并将该最优个体与种群中已知最优个体相比较，若当前分组下寻得的最优个体的目标函数值小于种群中已知最优个体的目标函数值，则用当前分组下寻得的最优个体替换种群中已知最优个体；否则，保持种群中已知最优个体不变。Compare the vehicle serial number k+1 with the number of vehicles M, if (k+1)≤M, set k=k+1 and return to step 3; otherwise, calculate the optimal individual found under the current group, and use the The optimal individual is compared with the known optimal individual in the population, and if the objective function value of the optimal individual found under the current group is smaller than the objective function value of the known optimal individual in the population, then the optimal individual found under the current group is used. The optimal individual replaces the known optimal individual in the population; otherwise, the known optimal individual in the population remains unchanged.

步骤7、上层终止判断Step 7. Upper layer termination judgment

将聚类调整次数I+1与上层调整次数U进行比较，若(I+1)≤U，则令I＝I+1并对输入的客户重新进行聚类，确定聚类中心(X_i,Y_i)，i＝1,2,…,M，并返回步骤2；否则，停止聚类，并返回记录的已知最优路径、路径长及迟到时间，该路径即为所选多时间窗车辆路径的最优路径。Compare the number of clustering adjustments I+1 with the number of upper-level adjustments U, if (I+1)≤U, set I=I+1 and re-cluster the input customers to determine the cluster center (X_i , Y_i ), i=1,2,...,M, and return to step 2; otherwise, stop clustering, and return the recorded known optimal path, path length and late time, the path is the selected multi-time window Optimal route for vehicle routing.

本发明的效果可通过下述实验进一步说明。The effects of the present invention can be further illustrated by the following experiments.

1、实验条件1. Experimental conditions

（1）实验数据来源(1) Source of experimental data

由于多时间窗车辆路径问题目前还没有标准的测试数据集，本发明的VRPMTW测试数据集是在单时间窗的VRPTW Solomon Benchmark测试数据集上，每个客户节点单独随机增加1个时间窗，生成规模如下：Since there is no standard test data set for the multi-time window vehicle routing problem, the VRPMTW test data set of the present invention is based on the VRPTW Solomon Benchmark test data set with a single time window, and each client node randomly adds a time window separately to generate The scale is as follows:

a.从VRPTW Solomon Benchmark的离散型数据集R1和聚焦型数据集C1中各取3组数据，为每组数据的每个客户节点增加1个与原来时间窗不重合的时间窗，同时保持其原有数据不变；a. Take 3 sets of data from VRPTW Solomon Benchmark's discrete data set R1 and focused data set C1 respectively, and add a time window that does not overlap with the original time window for each customer node of each set of data, while maintaining its The original data remains unchanged;

b.若原有时间窗长度超过100，则新产生时间窗长度为0到100之间的随机数；否则，新产生时间窗与原时间窗长度相等。b. If the length of the original time window exceeds 100, the length of the newly generated time window is a random number between 0 and 100; otherwise, the length of the newly generated time window is equal to the length of the original time window.

（2）实验参数(2) Experimental parameters

种群规模：p_s＝25           客户数目：n_c＝50或100Population size: p_s =25 Number of customers: n_c =50 or 100

聚类次数：U＝5             聚类调整次数：D＝5Clustering times: U=5 Clustering adjustment times: D=5

进化迭代次数：I_sub＝200    服务车辆数：5＜M＜15Number of evolution iterations: I_sub = 200 Number of service vehicles: 5<M<15

种群施力个体数：m从种群规模逐渐减小到1。The number of individuals exerting force in the population: m gradually decreases from the population size to 1.

2、实验内容2. Experimental content

（1）离散型多时间窗车辆调度问题的测试(1) Test of discrete multi-time window vehicle scheduling problem

分别用现有的万有引力算法GSA和本发明方法IGSA对离散型数据集R1进行测试，记录各自的最优路径长度、最差路径长度、路径长均值、路径长标准差和达到时间均值，结果如表1所示：Existing gravitational algorithm GSA and the method IGSA of the present invention are used to test the discrete data set R1 respectively, and record respective optimum path length, worst path length, path length mean value, path length standard deviation and arrival time mean value, the result is as follows Table 1 shows:

表1R1数据集测试结果Table 1R1 data set test results

（2）聚集型多时间窗车辆调度问题的测试(2) Test of aggregated multi-time window vehicle scheduling problem

分别用现有的万有引力算法GSA和本发明方法IGSA对聚集型数据集C1进行测试，记录各自的最优路径长度、最差路径长度、路径长均值、路径长标准差和迟到时间均值，结果如表2所示：Existing gravitational algorithm GSA and the method IGSA of the present invention are used to test aggregated data set C1 respectively, and record respective optimum path length, worst path length, path length mean value, path length standard deviation and late time mean value, the results are as follows Table 2 shows:

表2C1数据集测试结果Table 2C1 Data set test results

3、实验分析3. Experimental analysis

从表1可知，对客户规模为50的数据集R102、R103及客户规模为100的数据集R102、R103的求解结果表明，本发明方法求得的结果无迟到或迟到时间在10左右，即求得结果客户满意度较高。客户规模为50的数据集R101及客户规模为100的数据集R101的求解结果表明，由于时间窗较窄，求得结果的迟到时间有明显增加，均在100左右，即求得结果客户满意度有下降。但比较最优路径长度、最差路径长度、路径长度均值及路径长标准差这4项的数据，可得，本发明方法求得的结果都一致小于现有的万有引力算法求得的结果。因此，在对离散型的多时间窗车辆路径问题进行求解时，本发明方法比现有的万有引力算法性能更优。As can be seen from Table 1, the solution results of the data sets R102, R103 with a customer scale of 50 and the data sets R102 and R103 with a customer scale of 100 show that the results obtained by the method of the present invention are not late or the late time is around 10, that is, The result is high customer satisfaction. The solution results of the data set R101 with a customer scale of 50 and the data set R101 with a customer scale of 100 show that due to the narrow time window, the late arrival time of the obtained results has increased significantly, both of which are around 100, that is, the customer satisfaction of the obtained results There is a decline. But comparing the data of these 4 items of optimal path length, worst path length, path length mean value and path length standard deviation, it can be obtained that the results obtained by the method of the present invention are all consistently less than the results obtained by the existing gravitational algorithm. Therefore, when solving the discrete multi-time window vehicle routing problem, the performance of the method of the present invention is better than that of the existing universal gravity algorithm.

从表2可知，对客户规模为50的数据集C102、C103及客户规模为100的数据集C103的求解结果表明，本发明方法求得结果均无迟到，现有的万有引力算法求得的结果几乎无迟到，即求得结果客户满意度较高。对客户规模为50的数据集C101及客户规模为100的数据集C101、C102的求解结果表明，本发明方法和现有的万有引力算法求得的结果均出现了10⁺²数量级迟到，即客户满意度下降，因此对时间窗窄的数据集，要求的客户满意度更高的路径还需要调整加大扫描搜索的细度和迭代次数。比较最优路径长度、最差路径长度、路径长度均值及路径长度标准差这4项的数据，可得，本发明方法求得的结果都优于现有的万有引力算法。因此，在对聚集型的多时间窗车辆路径问题进行求解时，本发明方法比现有的万有引力算法性能更优。As can be seen from Table 2, the solution results for data sets C102 and C103 with a customer scale of 50 and data set C103 with a customer scale of 100 show that the results obtained by the method of the present invention are not late, and the results obtained by the existing gravitational algorithm are almost No lateness means high customer satisfaction. The results of solving the data set C101 with a customer scale of 50 and the data sets C101 and C102 with a customer scale of 100 show that the results obtained by the method of the present invention and the existing gravitational algorithm all appear to be 10⁺² orders of magnitude late, that is, customer satisfaction Therefore, for data sets with narrow time windows, the path that requires higher customer satisfaction needs to be adjusted to increase the fineness of scanning search and the number of iterations. Comparing the data of the four items of optimal path length, worst path length, path length mean value and path length standard deviation, it can be obtained that the results obtained by the method of the present invention are better than the existing gravitational algorithm. Therefore, when solving the aggregated multi-time window vehicle routing problem, the performance of the method of the present invention is better than that of the existing universal gravity algorithm.

综上，在对离散型和聚集型的多时间窗车辆路径问题进行求解时，本发明方法比现有的万有引力算法具有更强的搜索能力，有效性和实用性突出。To sum up, when solving discrete and aggregate multi-time window vehicle routing problems, the method of the present invention has stronger search ability than the existing universal gravity algorithm, and its effectiveness and practicability are outstanding.

Claims (3)

1. A multi-time window vehicle path selection method based on an improved universal gravitation algorithm comprises the following steps:

(1) inputting the distribution position condition of a customer, the number M of vehicles in a parking lot, the lower-layer iteration times D and the upper-layer adjustment times U, and carrying out rapid clustering processing on the customer, namely randomly selecting M points as a clustering center, and setting the clustering adjustment times I to be 1;

(2) clustering and grouping: numbering each vehicle, and enabling a serial number k of a first vehicle to be 1; dividing a customer distribution area into a plurality of vehicle areas by taking a vehicle yard as a center according to a clustering center, and determining a customer set served by each vehicle;

(3) initializing a population: recording a client set served by a kth vehicle as a group, and initializing the spatial position of each individual in the group at the moment when t is 0;

(4) optimizing by utilizing a universal gravitation algorithm: at the time T, selecting m individuals with the best fitness value in the population, applying force to other individuals in the population under the action of universal gravitation, then updating the individual speed and the individual position once to obtain the initial speed and position of the population at the time T +1, and enabling T to be T + 1;

(5) local search:

(5a) performing one-time cross operation on each individual in the updated population to obtain a temporary population, performing boundary constraint check on each individual in the temporary population, and replacing the initial individual with the individual with a target value superior to that of the initial individual if the target value of the individual exists in the population, otherwise, keeping the position of each individual in the population unchanged;

(5b) comparing T with the lower layer iteration times D, and if T is less than or equal to D, returning to the step (4); otherwise, recording the optimal path, the path length and the late time of the kth vehicle;

(6) and (3) lower layer termination judgment: comparing the vehicle serial number k +1 with the vehicle number M, and if (k +1) is less than or equal to M, returning to the step (3); otherwise, calculating the optimal individual searched under the current grouping, comparing the optimal individual with the known optimal individual in the population, and if the objective function value of the optimal individual searched under the current grouping is smaller than the objective function value of the known optimal individual in the population, replacing the known optimal individual in the population with the optimal individual searched under the current grouping; otherwise, keeping the known optimal individuals in the population unchanged;

(7) and (3) upper layer termination judgment: comparing the clustering adjustment times I +1 with the upper-layer adjustment times U, if the (I +1) is less than or equal to U, re-clustering the input clients, determining a clustering center, and returning to the step (2); otherwise, stopping clustering, and returning the recorded known optimal path, path length and late time, wherein the path is the optimal path of the selected multi-time window vehicle path.

2. The improved gravity algorithm-based multi-time window vehicle routing method according to claim 1, wherein the step (3) initializes the spatial position of each individual in the population according to the following formula:

<math> <mrow> <msubsup> <mi>p</mi> <mi>j</mi> <mi>i</mi> </msubsup> <mo>=</mo> <mrow> <mo>(</mo> <mi>b</mi> <mo>-</mo> <mi>a</mi> <mo>)</mo> </mrow> <mi>R</mi> <mo>+</mo> <mi>a</mi> <mo>,</mo> <msubsup> <mi>p</mi> <mi>j</mi> <mi>i</mi> </msubsup> <mo>&Element;</mo> <mo>[</mo> <mi>a</mi> <mo>,</mo> <mi>b</mi> <mo>]</mo> <mo>,</mo> <mi>i</mi> <mo>=</mo> <mn>1,2</mn> <mo>,</mo> <mo>.</mo> <mo>.</mo> <mo>.</mo> <mo>,</mo> <mi>L</mi> <mo>,</mo> <mi>j</mi> <mo>=</mo> <mn>1,2</mn> <mo>,</mo> <mo>.</mo> <mo>.</mo> <mo>.</mo> <mo>,</mo> <msub> <mi>N</mi> <mi>k</mi> </msub> </mrow></math>

wherein,

is the position coordinate of the ith individual in the jth dimension space, and R is [0,1 ]]Random number in, L is the size of the population, k represents the number of the vehicle, N_kIndicates the number of customers served by the k-th vehicle, [ a, b ]]Is a spatial range.

3. The improved gravity algorithm-based multi-time window vehicle routing method according to claim 1, wherein the updating of the individual speed and the individual position in step (4) is performed once according to the following formula:

<math> <mrow> <msubsup> <mi>v</mi> <mi>i</mi> <mi>d</mi> </msubsup> <mrow> <mo>(</mo> <mi>t</mi> <mo>+</mo> <mn>1</mn> <mo>)</mo> </mrow> <mo>=</mo> <mi>R</mi> <mo>&times;</mo> <msubsup> <mi>v</mi> <mi>i</mi> <mi>d</mi> </msubsup> <mrow> <mo>(</mo> <mi>t</mi> <mo>)</mo> </mrow> <mo>+</mo> <msubsup> <mi>a</mi> <mi>i</mi> <mi>d</mi> </msubsup> <mrow> <mo>(</mo> <mi>t</mi> <mo>)</mo> </mrow> <mo>,</mo> <msubsup> <mi>p</mi> <mi>i</mi> <mi>d</mi> </msubsup> <mrow> <mo>(</mo> <mi>t</mi> <mo>+</mo> <mn>1</mn> <mo>)</mo> </mrow> <mo>=</mo> <msubsup> <mi>p</mi> <mi>i</mi> <mi>d</mi> </msubsup> <mrow> <mo>(</mo> <mi>t</mi> <mo>)</mo> </mrow> <mo>+</mo> <msubsup> <mi>v</mi> <mi>i</mi> <mi>d</mi> </msubsup> <mrow> <mo>(</mo> <mi>t</mi> <mo>+</mo> <mn>1</mn> <mo>)</mo> </mrow> </mrow></math>

wherein,

(t) is the speed of the individual i in the d-dimensional space at time t,(t) acceleration of the individual i in the d-dimensional space at time t,

(t) the position of the individual i in the d-dimensional space at time t, R is [0,1 ]]The random number in (c).

Images