Abstract
In the context of sharing economy, logistics companies have begun to adopt a new collection and distribution model based on external vehicles, and external vehicles are used to provide customers with collection/distribution services. A kind of Multi-Vehicle Split Pickup and Delivery Problem with Time Windows (MVSPDPTW) is studied in the paper. The minimum total length of the vehicle’s travel path is the objective function, and a mixed-integer linear programming model is established. A high-efficiency Tabu Simulated Annealing (TSA) algorithm is proposed. Two new neighborhood search operators are designed in the algorithm, they are used to repair the violation of capacity constraints and the operation of car replacement. In the method, the neighborhood search range is expanded and the algorithm is avoided from falling into a local optimum. In addition, a taboo mechanism and a penalty mechanism for violation of constraints are added to the algorithm, the effective tailoring of the search space is realized and the algorithm’s global optimization ability is improved. Finally, a large number of simulation experiments were performed on the algorithm based on the Solomon data set and the constructed simulation data set, and the effectiveness of the algorithm is verified in the experiments.
Introduction
The development of network economy is inseparable from logistics technology, and logistics development is inseparable from timeliness and cost saving. In the sharing economy mode, the idea of this paper is to study the optimal allocation of vehicle material receiving and dispatching under this environment, the utilization rate of social resources is improved.
Traditional logistics companies use their own vehicles to provide customers with collection/distribution services. As self-owned vehicles have disadvantages such as high cost of use and difficult management, in the context of the sharing economy, logistics companies have begun to adopt a new collection and distribution model based on external vehicles, and external vehicles are used to provide customers with collection/distribution services. The external vehicles are mainly social citizen idle vehicles. The new collection and distribution model can reduce the use of company-owned vehicles, reduce logistics costs, and it is of great significance to improving the utilization rate of social citizen idle vehicles and environmental protection.
Based on the new pickup and delivery model, a Multi-Vehicle Split Pickup and Delivery Problem with Time Windows (MVSPDPTW) is studied in this paper, social citizen idle vehicles are used as the resources for the pickup and delivery service. For collection/distribution tasks, social citizen idle vehicles have their own starting point, destination, earliest departure time and latest return time. Vehicles must depart from the starting point after the earliest departure time to provide customers with collection/distribution services, and return to the end at the latest before the return time. In the process of serving customers, vehicles need to provide services within the customer’s time window and allow demand splitting. The same customer demand can be fulfilled by multiple vehicles in order to make full use of the vehicle’s loading capacity.
Vehicle Routing Problems (VRP) have been a hot issue in the field of logistics research since it was proposed in 1959 [1, 2]. At present, in the research of VRP, it is generally restricted that the customer’s needs can only be fulfilled by one car, it belongs to the VRP whose demand cannot be separated. In 1989, Split Delivery Vehicle Routing Problems (SDVRP) with split demand was studied for the first time [3, 4]. Reasonable splitting of customer needs was beneficial to increase the vehicle loading rate and save cost. In 1994, these were proved that SDVRP is an NP-hard problem [5]. Many scholars have conducted research on SDVRP from different angles. Since customers usually require services within a given time, a time window-constrained SDVRP has been studied. The time limit is introduced to SDVRP for the first time [6], a time-window-constrained SDVRP model is constructed, and a two-stage heuristic algorithm is proposed for solving. A tabu search algorithm is proposed to solve the SDVRP with time window constraints [7], and a variety of new neighborhood search operators are designed. SDVRP with soft time window constraints is studied [8–11]. The soft time window constraints allow customers to be served outside a given time window by adding a penalty function. In SDVRP, only the customer’s distribution needs are considered. With the development of reverse logistics, customers often have collection/distribution requirements at the same time. Therefore, the SDVRP of a collection demand has been studied, and the demand can be split into the collection and distribution problem [12, 13], this reveals the benefits of splitting requirements. For a given set of collection and distribution locations, when the customer demand is slightly more than half of the vehicle’s load, it proves to bring the greatest benefit. An improved ant colony algorithm is proposed to solve the demand-separable collection and distribution [14, 15], and the experimental results show that the algorithm can improve the operation efficiency of the distribution center.
In the past, the research on SDVRP was basically based on a single vehicle type, but in actual distribution, vehicles often have multiple vehicle types. The SDVRP of multiple vehicle types is studied, and the loading capacity of the vehicle and the maximum mileage are different, and a simulated annealing algorithm is proposed to solve it [16]. Open vehicle routing problem (OVRP) is another type of vehicle routing problem. The main difference between it and the basic vehicle routing problem (VRP) is that the vehicle is not required to return to the original starting point after completing the pickup and delivery task, or if it is required to return to the original starting point, it will return along the original route. Therefore, the driving route of vehicles is open rather than closed. The multi-model open vehicle routing problem with splittable requirements is studied, the vehicle starts from different distribution centers to serve customers, and ends the service after serving the last customer without returning to the distribution center, a tabu search algorithm is proposed to solve the problem [17]. The SDVRP of single vehicle type and multiple vehicle types are studied separately [18], the experimental results show that vehicles with different loading capacities are more in line with actual needs, and it can effectively increase the loading rate of vehicles.
Most of the above-mentioned documents separately study the demand for collection, distribution, multiple models, and customer time windows. However, with the rapid development of logistics technology, user needs are more complex, and this single-mode research has been difficult to meet actual application needs. Therefore, the above-mentioned multiple factors are comprehensively considerd in this article, and a multi-vehicle demand with a time window is researched, the demand can be split to collect and distribute. To solve this problem, the minimum total length of the vehicle path to perform the task is taken as the objective function, and a mixed integer linear programming model is established, an efficient Tabu Simulated Annealing (TSA) algorithm is proposed to solve the problem.
Mathematical model
Briefly describe the process
In MVSPDPTW, there is a given distribution center, a group of customers with different geographical locations and a group of vehicles with multiple models. The goal of MVSPDPTW is to select the right vehicle to perform all customer tasks, so that the total length of the vehicle’s travel path is minimized.
0 means the distribution center, T ={ 1, 2, …, n } means the collection of customer points, each customer point corresponds to a collection/distribution task, the customer point corresponding to the collection task is called the pickup point, the customer point corresponding to the delivery task is called the delivery point. For customer point i ∈ T, [e i , l i ] represents the i-th time window, e i represents the earliest service start time, and l i represents the latest service start time. q i represents the quantity of goods i need to collect/deliver. Each customer point can be visited by multiple vehicles, but it is only allowed to be visited by the same vehicle at most once. K = { 1, 2, …, m } represents the set of vehicles. For k ∈ K, s k , n k , Q k , o k , and f k represent the starting point, ending point, maximum capacity, earliest departure time, and latest return time of vehicle k, respectively. As shown in Fig. 1, the vehicle needs to depart from the starting point after the earliest departure time, collect the goods at the pickup point and send them back to the distribution center, load the goods at the distribution center and deliver them to the delivery point, and finally return to the destination before the latest return time. The vehicle is allowed to visit the distribution center multiple times during the mission.

Schematic diagram of vehicle execution tasks.
Since the vehicle may visit the distribution center many times, if a variable is used in the model to distinguish how many times the vehicle visits the distribution center, the difficulty of modeling and solving will increase significantly. In order to reduce the difficulty, the following processing methods are adopted: T ={ 1, 2, …, n } represents the collection of customer points, n is the actual number of customer points, for the pickup point numbered i, a numbered n + i virtual delivery point is constructed to correspond to it; for the delivery point numbered j, it is renumbered as n + j, and a virtual pickup point numbered j is constructed to correspond to it. After the above processing, T = T P ∪ T D , T P ={ 1, 2, …, n } represents the set of pickup points, T D ={ n + 1, n + 2, …, 2n } represents the set of delivery points, and the vehicle needs to collect the goods from point i (i ∈ T P ) and deliver them to point n + j (n + j ∈ T D ). In this processing method, the distribution center is replaced with n virtual customer points, all virtual customer points are located at the same location as the distribution center and there is no service time window restriction; the amount of goods collected/delivered by each virtual customer point corresponds to the equal actual customer in the points. As shown in Fig. 2, the distribution center in Fig. 1 is virtualized and replaced with virtual customer points in the box, where i* = n + j.

Virtualization schematic.
After the above processing, MVSPDPTW is defined on the upper graph of the fully directed graph G = (V, A), V ={ S ∪ N ∪ T ∪ 0 } represents the vertex set, A ={ (i, j) : i, j ∈ V, i ≠ j } represents the arc set.
S: The starting point set of the vehicle, S ={ s1, s2, …, s m };
N: the set of end points of the vehicle, N ={ n1, n2, …, n m };
T: customer point collection, T ={ 1, 2, …, 2n };
d ij : the distance between point i and point j, i, j ∈ V;
t ij : the time between point i and point j, i, j ∈ V;
b ik : the time when the vehicle k arrives at the point i, k ∈ K, i ∈ T ∪ 0;
a ik : the start time of vehicle k at point i, k ∈ K, i ∈ T ∪ 0;
r ik : the service time of vehicle k at point i, k ∈ K, i ∈ T ∪ 0;
q ik : the amount of cargo loaded/unloaded by vehicle k at point i, k ∈ K, i ∈ T ∪ 0;
w ik : Cargo capacity of vehicle k after point i, k ∈ K, i ∈ T ∪ 0;
x ijk : Indicates whether the vehicle k is directly from point i to point j, if yes, it is 1, otherwise it is 0;
M: Represents a very large integer.
Linear model
Wherein, formula (1) is the objective function, which minimizes the total length of the mission vehicle’s travel path; formulas (2) and (3) indicate that the vehicle will perform the task from the starting point and finally return to the end; formulas (4) and (5) indicate every customer point can be accessed by multiple vehicles, but it is only allowed to be accessed by the same vehicle at most once; formulas (6), (7), (8) indicate that the vehicle needs to collect goods from the pickup point for delivery to the corresponding delivery point; formula (9) means that the cargo carrying capacity of the vehicle from the starting point and returning to the destination is 0; formula (10) expresses the capacity constraint; formulas (11) and (12) express the change of the vehicle carrying cargo capacity at the customer point; formula (13) means the sum of the divided cargo volume of each customer point is equal to the total cargo volume of the customer point; formulas (14) and (15) indicate that the vehicle needs to depart from the starting point after the earliest departure time to provide services to customers, and return before the latest return time end point; formulas (16) and (17) express the relationship between the time when the vehicle arrives at the customer point and the time when it starts service at the customer point; formula (18) expresses that the vehicle can only provide service within the customer’s time window.
In order to solve MVSPDPTW, a tabu simulated annealing algorithm is proposed, and two new neighborhood search operators are designed, which are used to repair the violation of vehicle capacity constraints and the operation of changing vehicles. A taboo mechanism and a penalty mechanism are added for violation of constraints in the algorithm, the effective tailoring of the search space is realized, and the algorithm’s global optimization capability is improved. First, the dynamic greedy algorithm is used to generate the initial feasible solution, and the initial feasible solution is used as the input of the tabu simulated annealing algorithm, and finally the global best solution is output.
Initial solution generation
The simulated annealing algorithm is essentially an iterative solution algorithm. The initial solution determines the starting point of the iteration. Therefore, a good initial solution helps the algorithm find a better final solution. In the initial solution generation method in this article, a simple and effective dynamic greedy algorithm is used to generate the initial solution [5].
Formulas (20) and (21) respectively represent the selection criteria of vehicles and customer points, and the selection of vehicles is based on the principle of maximizing service time. The selection of customer points is based on the principle of minimizing the total travel time and waiting time between the last visit object i in the current planning path and the customer point j to be selected. Formula (21) is used to select service customer points for the current vehicle in turn, and to check whether the selected customer points meet the insertion constraints, including time window constraints, capacity constraints, customer point visit times constraints, and demand constraints. Demand constraints refer to all delivery of the goods which come from the distribution center, and all collected goods are sent back to the distribution center. If the constraints are met, the customer points are inserted into the current vehicle planning path, and this process is repeated until all customer needs are met.
The specific process of the dynamic greedy algorithm is as follows:
In the iterative process, accepting infeasible solutions can help increase the search range of the neighborhood and improve the ability of the algorithm to jump out of the local optimum, and it helps the algorithm to obtain a better feasible solution in the iterative step [19], so a violation constraint punishment mechanism is designed in this paper. The overall constraint penalty function is designed in formula (22), B represents the fitness value of the path planning solution, β represents the parameter penalty weight, η (i, k) , φ (i, k) , λ (i, k) represent the capacity constraint evaluation value, time window constraint evaluation value and demand constraint evaluation value of point i in R k , respectively. Their calculation methods are given by the formula (23), (24), (25). R k ={ o k , …, n k } represents the execution task path of vehicle k, which is arranged according to the order of vehicle visits.
Tabu search is a meta-heuristic algorithm that simulates the human thought process. Its basic idea is to avoid repeated searches in the searched space and force the algorithm to explore other new solutions.
The current solution of each iteration is taken as the taboo object and stored in the tabu list [20]. In the subsequent iterations of l (the length of the tabu list), the object is set to the taboo state, which is different from the side movement taboo in taboo structures [21, 22]. The length of the tabu list refers to the number of objects which are stored in the tabu list. When the number of the stored objects reaches the upper limit and another object is added, the first in, first out principle is adopted to remove the previous objects from the tabu list and then add them.
Multi neighborhood structure
The multi-neighbor structure can increase the diversity of solutions and improve the possibility of the algorithm to find the best global solution. In this paper, five kinds of neighborhood search operators are used, the migration operator, the migration division operator, and the exchange operator are derived from the classic neighborhood search operator [23]. In the iterative process of this article, the repair operator and the car exchange operator are designed to accept the infeasible solution and the characteristics of multiple models, and they are used to repair the violation of capacity constraints and the operation of changing vehicles respectively.
In order to better explain the operator operation, the following symbols are defined. R
k
(k ∈ K) and R
h
(h ∈ K) represent the execution task path of vehicle k and vehicle h, respectively. For client points

Schematic diagram of the task set.

Schematic diagram of migration operator operation.

Schematic diagram of migration and division operator operation.

Schematic diagram of exchange operator operation.

Schematic diagram of repair operator operation.

Schematic diagram of car change operator operation.
Since a single customer point is only allowed to be accessed by the same vehicle at most once, after the operator operation, if the situation shown in Fig. 9 occurs, R k contains two customer points i, then one i is selected to be deleted from R k based on the principle of minimizing B (R k ), and to combine the quantity of goods which are loaded/unloaded by two services i of k.

Schematic diagram of merged customer points.
The dynamic greedy algorithm is used to generate the initial feasible solution, and it is set to the current solution S now of the 0-th generation and the global best solution S best , and then the iteration starts. In each iteration, the design criteria is used to update S now and S best . Finally, when the current temperature is less than the set minimum temperature, the algorithm terminates and outputs S best . The specific algorithm flow is as follows:
If e(Scan-Snow)/T > Random (0, 1), S can will become the new S now , Random(0,1) is used to generate random data uniformly distributed on (0,1), otherwise keep S now .
Simulation and performance analysis
In this paper, the VRP (Vehicle Routing Problem with Time Windows, VRPTW) data set with time window constraints [24], the VPRSD (Split Delivery Vehicle Routing Problem with Time Windows, SDVRPTW) data set with time window constraints [25] and the constructed MVSPDPTW simulation data set are used to test TSA algorithm, the experimental results were compared and analyzed. Each example was run 20 times. All experiments were performed on a desktop computer with i5-8500 CPU 3.0 GHz, 16 GB memory, and myeclipse 2017 environment.
Parameter setting
Setting reasonable parameters can effectively improve the performance of the algorithm [26]. First, the parameter settings in references [27, 28] are used to determine the test interval of each parameter in this article, and then different scale examples are selected to perform multiple tests by changing the parameter combination. Finally, the best combination of parameters is determined.
The value range and step length of each parameter are shown in Table 1. In order to simplify the calculation complexity, the controlled variable method is used to determine the optimal value of each parameter in turn. The final result of TSA parameter setting is as follows: initial temperature T0 = 100, termination temperature T end = 0.1, iteration length L = 30, cooling rate R T = 0.93, penalty weight β= 4, number of neighborhood solutions N u = 2n + 150 (n is the number of customers), the length of the tabu list is l = 8.
Algorithm parameter combination
Algorithm parameter combination
As there is currently no recognized data set available in the field of MVSPDPTW, there is no way to directly test the algorithm in this article. MVSPDPTW is an extension of the classic VRPTW and SDVRPTW. Many of the constraints in the problem are similar. Therefore, if the algorithm in this article can solve VRPTW and SDVRPTW well, then there is reason to believe that the algorithm in this paper can also effectively solve MVSPDPTW.
Numerical experiment of VRPTW dataset
From the VRPTW data set, 27 VRPTW examples are selected for testing [24], and the vehicle capacity is all 200 in the VRPTW examples. Table 2 records the comparison between R best and the known best solution K best of the TSA algorithm for each example. Gap represents the gap between R best and K best , and Gap = (R best - K best )/K best . K best data comes from http://web.cba.neu.edu/ msolomon/problems.htm. The design highlights several factors that affect the behavior of routing and scheduling algorithms. They are: geographical data; the number of customers serviced by a vehicle; percent of time-constrained customers; and tightness and positioning of the time windows. The geographical data are randomly generated in problem, Problem sets have a short scheduling horizon and allow only a few customers per route (approximately 5 to 10). In contrast, there are a long scheduling horizon permitting many customers (more than 30) to be serviced by the same vehicle. The customer coordinates are identical for all problems within one type C. The problems differ with respect to the width of the time windows. Some have very tight time windows, while others have time windows which are hardly constraining. In terms of time window density, that is, the percentage of customers with time windows, I created problems with 25, 50, 75 and 100 % time windows.
The TSA algorithm solution result is compared with the known best solution
The TSA algorithm solution result is compared with the known best solution
According to the data in Table 2, the TSA algorithm has updated 4 VRPTW examples (example number C101, number of customers 25; example number C102, number of customers 50; example number C107, number of customers 50; example number C109, number of customers 50) in the best known solution. In the remaining 23 VRPTW calculation examples, although the TSA algorithm’s solution results are inferior to the known best solution, the gap between them is kept within 1%. The experimental results show that the TSA algorithm is effective in solving VRPTW.
From the SDVRPTW data set [25], 12 SDVRPTW examples were selected for testing. In the SDVRPTW example, the vehicle capacity is all 100. Table 3 records the comparison between TSA algorithm and BPC algorithm (Branch and Price and Cut algorithm) for solving the 12 SDVRPTW calculation examples. Gap represents the gap between the TSA algorithm to solve each example result (R best ) and the BPC algorithm to solve each example result (C best ), Gap = (R best -C best C best )/ C best .
Comparison of TSA algorithm and BPC algorithm solution results
Comparison of TSA algorithm and BPC algorithm solution results
From the data in Table 3, it can be seen that the solution results of the TSA algorithm are the same as the BPC algorithmin the small-scale case with 25 customers. With the increase in the number of customers, although the solution results of the TSA algorithm are inferior to the BPC algorithm, the gap between them remains within 1%. In terms of solution time, the BPC algorithm solution time is less than the TSA algorithm in small-scale cases, but as the number of customers increases, the BPC algorithm solution time far exceeds the TSA algorithm. In the solution time for the same number of customers and different cases, the TSA algorithm remains basically unchanged, while the BPC algorithm is significantly different. The experimental results fully demonstrate that the TSA algorithm is close to the BPC algorithm in solving SDVRPTW, but the TSA algorithm takes less time.
12 MVSPSPTW simulation examples were constructed based on 12 SDVRPTW examples. Without changing the distribution center location, customer location, customer cargo volume, and customer time window in the SDVRPTW example, only half of the customers were randomly selected to change their distribution tasks to collection tasks. The starting point and ending point of the vehicle are no longer set as the distribution center, but they are randomly generated within the 100×100 square area that the customer belongs to. The starting service time of a vehicle is a randomly generated integer between 0 and 800, and the end service time is equal to a randomly generated integer between 500 and 900 plus the starting service time of the vehicle; the maximum vehicle capacity is a randomly generated integer between 70 and 120. In addition, when the number of customers is 25, 15 vehicles are randomly generated, and when the number of customers is 50, 30 vehicles are randomly generated to ensure sufficient vehicle resources.
In Table 4, the TSA algorithm, the simulation results’ comparison of the traditional simulated annealing (SA) algorithm, and particle swarm optimization (PSO) algorithm to solve the MVSPDPTW are recorded. The design and processing method of the PSO algorithm comes from the reference [29, 30]. The parameters of the traditional SA algorithm and the PSO algorithm are determined according to the parameter setting method in Table 1. Gap1 and Gap2 respectively represent the gap between TSA algorithm solving each example result (R best ) and SA algorithm solving each example result (T best ), PSO algorithm solving each example result (P best ), Gap1 = (T best - R best )/R best , Gap2 = (P best - R best )/R best .
Comparison of TSA, SA, PSO algorithm solution results
Comparison of TSA, SA, PSO algorithm solution results
According to the data in Table 4, the TSA algorithm in this paper has better solution performance than the traditional SA algorithm and PSO algorithm. The main factors are: (i) The TSA algorithm has the best solution result, which is up to 2.40% higher than the traditional SA algorithm, and there is an average increase of 1.14%. Compared to the PSO algorithm, the highest increase is 5.19%, and an average increase is 3.52%. (ii) In terms of solving time, the SA algorithm is the longest, the TSA algorithm is second, and the PSO algorithm is the shortest. However, the overall time-consuming gap between the TSA algorithm and the PSO algorithm is not large.
In addition, the required time for the TSA algorithm to solve the same number of customers is basically the same, and it has good stability. The main reasons for the better solution performance of the TSA algorithm are: (i) Due to the characteristics of the multi-models in this paper, the result of the car selection will have a vital impact on the pros and cons of the final solution. The car change operator is helpful to find the most suitable vehicle to perform the task. In addition, a variety of operator coordination methods are beneficial to expand the neighborhood search range and improve the algorithm’s global optimization capability. (ii) A taboo mechanism and a penalty mechanism for violation of constraints are added to the TSA algorithm, the effective tailoring is achieved in the search space, avoiding the algorithm from falling into a local optimum and saving solution time. In contrast, the traditional SA algorithm and the PSO algorithm are easy to fall into a local optimum.
In this paper, under the new mode of solicitation and distribution, we study a kind of multi-vehicle demand with a time window that can be divided into solicitation and distribution. To solve this problem, a mixed-integer linear programming model is established. In the process of constructing the model, a corresponding virtual collection/distribution point was constructed for each customer point, the increase of some variables is reduced in the model, and the big M method was used to linearize the nonlinear constraints. An efficient TSA algorithm is designed. In addition to the traditional SA algorithm framework, the TSA algorithm incorporates a taboo mechanism and a penalty mechanism for violation of constraints, which realizes the effective tailoring of the search space and improves the algorithm’s global optimization capability. The TSA algorithm updates 4 known best solutions in the 27 VRPTW examples. In the solution results of the 12 SDVRPTW examples, the difference between the TSA algorithm and the BPC algorithm is not more than 1%, and the solution time is shorter. As far as the 12 MVSPDPTW simulation examples are concerned, the TSA solution results are not inferior to the traditional SA algorithm and PSO algorithm, and the maximum increase is 2.40% and 5.19%, respectively. The experimental results fully prove the feasibility and effectiveness of the TSA algorithm in this paper to solve MVSPDPTW.
There are still many problems that need to be solved urgently in the research content of this paper. Future research work can be carried out from the following two aspects: (i) Consider the demand to be split according to orders instead of arbitrarily continuous split; (ii) Consider the use of a joint delivery method between external vehicles and self-owned vehicles.
Footnotes
Acknowledgments
This work was sponsored by Scientific Research Project (NO.20C1085) of Hunan Provincial Education Department, China.
