Heterogeneous-vehicle distribution logistics planning for assembly line station materials with multiple time windows and multiple visits

Aiming at distribution logistics planning in green manufacturing, heterogeneous-vehicle vehicle routing problems are identified for the first time with multiple time windows that meet load constraints, arrival time window constraints, material demand, etc. This problem is expressed by a mathematical model with the characteristics of the vehicle routing problem with split deliveries by order. A hybrid ant colony optimization algorithm based on tabu search is designed to solve the problem. The search time is reduced by a peripheral search strategy and an improved probability transfer rule. Parameter adaptive design is used to avoid premature convergence, and the local search is enhanced through a variety of neighborhood structures. Based on the problem that the time window cannot be violated, the time relaxation rule is designed to update the minimum wait time. The algorithm has the best performance that meets the constraints by comparing with other methods.


Introduction
The high efficiency of assembly line operation and the low efficiency of material distribution have become the main issues in intelligent manufacturing. Distribution strategies and disordered management are the main causes of low distribution efficiency. Traditional material distribution has some shortcomings, such as poor continuity, poor distribution efficiency and poor delivery timeliness. Smooth and efficient material distribution is a prerequisite to ensure the smooth and continuous production of an assembly line. Material distribution route planning is an important application of the vehicle routing problem (VRP). The problem is more complicated due to the complex constraints and many distribution tasks. Vehicle route optimization reasonably reduces the number of vehicles and the total running time, which is in line with the concept of green logistics. Route planning is the key to efficient material distribution. Considering the limitations of actual physical factors inside workshops, Ho and Liao (2009) conducted research on the dynamic scheduling of vehicles to avoid collisions. Choi and Lee (2002) proposed a dynamic material distribution method based on static material distribution to avoid material shortages or accumulation. Jin and Zhang (2016) studied route planning and vehicle scheduling in manufacturing systems based on meta-heuristics. Umar et al. (2015) proposed the integrated dynamic scheduling and routing of job and automatedguided vehicles in a flexible manufacturing systems environment. They also optimized makespan and vehicle travel time to generate suitable material distribution routes.

Warehouse
Assembly line There is no research on the material distribution problem of multiple visits to each station in different time windows for heterogeneous vehicles based on the assembly line constraints. Materials are distributed in batches and small amounts at high frequencies because the line-side can store only small amounts of materials. Materials are assembled and consumed in accordance with the normal rhythm, so there are different time window requirements for each delivery. Fig.1 shows the schematic diagram of the material distribution system in the assembly workshop. To the best of our knowledge, single-depot and heterogeneous-vehicle VRPs are identified for the first time with multiple time windows (SHVRPMTW) that meet load constraints, arrival time window constraints, material demand. This problem is expressed by a mathematical model with the characteristics of the VRP with split deliveries by order (VRPSDO) (Xia & Fu 2018a). The model is solved linearly to verify the correctness of the model. This problem has many different types of constraints, each of which can be challenging to consider, even in isolation. A hybrid ant colony optimization algorithm based on tabu search (HACO-TS) is designed to solve large-scale problems. The search time is reduced by a peripheral search strategy and an improved probability transfer rule. Parameter adaptive design is used to avoid premature convergence, and the local search optimization ability is enhanced through a variety of neighborhood structures and peripheral search strategies. In the post processing stage, variable neighborhood search (VNS) is used to further optimize the results. Finally, the forward relaxation rule is designed to update the minimum wait time for each route in the solutions. Smooth material flow is an important condition for enterprises to realize just-in-time production. Based on the problem that the time window cannot be violated; the time relaxation rule is designed to update the minimum wait time.
The rest of this article is arranged as follows: Section 2 reviews the related literature. Section 3 describes the problems and establishes the mathematical model. Section 4 describes the improved algorithm in detail. In Section 5, the mathematical model is solved to verify the correctness. The results of the improved algorithm are compared. Section 6 analyzes and summarizes the research.

Related works
The VRP has attracted wide attention since it was first proposed by Dantzig and Ramser (Dantzig & Ramser, 1959). Many achievements have been made in addressing the VRP and its expansion problems, such as VRPs with time windows (VRPTW) (Paradiso et al., 2020), VRPs with split deliveries (SDVRP) (Archetti & Speranza, 2012), the fleet size and mixed VRP (FSMVRP) (Simi et al. 2015), VRPs with simultaneous delivery and pickup (VRPSDP) (Simsir & Ekmekci, 2019). Inspired by the VRP, some scholars have turned their attention to material distribution route planning in manufacturing systems, which is an important application of the VRP. The goal is to arrange the distribution route reasonably to complete tasks at the minimum cost and traveling distance under the known material plan. This problem is time sensitive and depends on the material plan. Dondo and Cerda (2007) established a mixed-integer programming model of the multidepot heterogeneous fleet VRP with time windows and proposed a three-stage heuristic algorithm to solve the problem. Fallahi et al. (2008) established an integer programming model for material distribution route optimization based on the actual situation of an automobile manufacturing company with the goal of minimizing the total traveling distance. Chiang et al. (2014) established a mathematical model to minimize the number of vehicles and the total traveling distance. They adopted a knowledge-based evolutionary algorithm to solve the problem. Fazlolahtabar et al. (2013) considered the triple criteria to determine the AGV optimal route. They proposed intelligent agents to apply mathematical programming methods to solve conflict deadlocks in the decision-making process.
This problem is a VRP with multiple time windows (VRPMTW). VRPMTW is not a new problem, although it has not been studied extensively. The study begins with Favaretto et al. (2007) first proposing the VRPMTW with multiple visits. They allowed fleets to visit the same customer multiple times within multiple time windows. An improved ACO was proposed to solve the problem. Belhaiza et al. (2014) studied a new hybrid VNS heuristic for the VRPMTW with the goal of the minimum waiting time and delay time. They proposed a recursive method to calculate the optimal allocation of waiting time along the route. Belhaiza et al. (2017) designed a new optimization algorithm based on hybrid genetic variable neighborhood search to optimize the VRPMTW. Belhaiza (2018) used a hybrid VNS algorithm for multicriterion optimization within the framework of game theory. A Pareto nondominant solution was selected from the search space of feasible solutions that satisfy a set of Nash equilibrium conditions. Beheshti et al. (2015) studied a multiobjective VRP with multiple prioritized time windows. They proposed a cooperative coevolutionary quantum-genetic algorithm to solve this problem. Larsen and Pacino (2019) designed an adaptive large neighborhood search framework to solve the VRPMTW. They proposed a novel solution that represents an incremental assessment to obtain the insertion cost. Bogue et al. (2020) proposed column generation and a post optimization heuristic algorithm that provided upper and lower bounds for the optimal solution cost of the VRPMTW. There is another type of VRPMTW in which a single time window must be selected from multiple time windows for each vehicle. Multiple time windows add complexity to the problem. To effectively determine the optimal start time for vehicle service, Hoogeboom et al. (2020) presented an exact polynomial time algorithm. Heterogeneous vehicle routing is another important application of VRP, which usually considers a limited or an unlimited fleet of capacitated vehicles, including the fleet size and mix VRP (FSMVRP) introduced by Golden et al. (1984) and the heterogeneous fixed fleet VRP (HFFVRP) introduced by Taillard and E(1999). There are variations in time window constraints on this problem. The problem has been studied for a long time, and fruitful results have been produced. To minimize the distance traveled and the number of vehicles, Ho et al. (2008) developed two hybrid genetic algorithms to solve the VRP with multiple depots. Salhi Said et al. (2013) proposed a heuristic algorithm based on set segmentation to solve the FSMVRP with return pickup. Matei et al. (2015) studied HFFVRP with minimized total cost. A genetic algorithm based on migration was combined with a local search program to solve the problem. Goeke and Schneider (2015) studied the electric VRPTW and mixed fleets in combination with a realistic energy consumption model that incorporates speed, gradient and cargo load distribution. They developed an adaptive large neighborhood search algorithm to solve the problem. Based on the constraint that depots must provide services to customers in the fuzzy time windows, Adelzadeh M et al. (2014) studied the multi-depot vehicle routing problem with fuzzy time windows and heterogeneous vehicles. Guezouli and Abdelhamid (2017) proposed a decision support system for multidepot FSM-VRPTW and provided a solution for multicriteria based on a genetic algorithm. Fachini and Armentano (2020) provided the first exact algorithms for the standard HFFVRPTW, and their effectiveness has been analyzed in detail. The research results of VRP and its extensions are rich, especially the algorithm design. However, in the assembly shop, the complexity of materials and frequent distribution increase the scale of the problem. The distribution planning problem with multiple time windows and multiple visits is the focus of the research. Moreover, the combination of VRPMTW and heterogeneous fleet constraints has not previously been studied, but it has important guiding significance for logistics distribution.

Problem description
With assembly workshops as the background, material distribution is pulled by the actual line sequence of the master production schedule (MPS) to carry out orderly operations. Material requirement planning (MRP) can be obtained by integrating MPS, manufacturing Bill of Material (BOM), product cycle time. MRP determines the material supply plan, which provides the exact time windows and quantity of the demands for each station. Then the material removal list is transferred to the warehouse for removal. The cycle time and work process of the production line have been completed in the planning layer. According to the logic of delivery time and quantity, the planning layer summarizes the demands to generate the distribution on-line planning. Therefore, when analyzing logistics optimization from the perspective of the logistics operation layer, the key is on the distribution logistics planning. Orders and the units of material demands are input data to be determined. Under a reasonable production plan, efficient vehicle distribution routes are designed to ensure that materials are delivered to stations in time based on workshop constraints. Delays in distribution are likely to result in production stoppages, which means high real-time requirements for materials. Only a small amount of material is temporarily stored in the station line-side to support production in a short time. Materials are delivered in a single package that makes maximum use of the line-side capacity. To meet different distribution needs, there are generally many types of vehicles in the workshop. Due to the limited workshop area, the vehicles are not subject to maximum time or maximum distance constraints. The problem satisfies the following assumptions: (1) The vehicle will not stop due to failure or traffic jams during the transportation.
(2) A certain amount of material has already been stored in the line-side before starting production. Material is consumed evenly during production. To ensure the continuity of production, the logistics system distributes materials to every station at intervals of a certain period of time. Therefore, a station receives materials in multiple time windows. Fig. 2 shows the curve of the line-side inventory level of a certain material over time, which is also the analysis of multiple time windows. For the stations with multiple time windows, each time window is associated with a virtual station except for the first time. Materials need to be transported to the station on the road in advance, so the demand signal is sent at the inventory (time point ). It can be seen from Fig.2, the station can accept distribution services only when the line-side inventory are consumed to less than or equal to the inventory (time point or ). or is the latest time for material distribution when material in the line-side is just completely consumed. After unloading service, the existing stock cannot exceed the max capacity , as , , and shown in Fig. 2. The arrival time earlier than needs to wait for service, because it will lead to the accumulation of materials in the line-side. The formulas for multiple demand time windows of the stations are designed as follows, where, is the average consumption rate of material i. is the materials contained in a single package. is the line-side inventory at time t.
is the safe inventory of the station. The unloading time is not taken into account in Fig. 2, meaning the inventory could be replenished instantaneously when materials arrive at the station.
In [ , ], the materials can be unloaded as . The unloading service is refused if the arrival time exceeds because there is a shortage of materials, even stop production. Assuming stable production, the consumption rate is constant. Regardless of the unloading time point in , , , and , the third production cycle is on R3 or its extension line, these moments do not affect the next demand window.

Mathematical model
Because of the complexity of the constraint of multiple time windows and multiple visits, a virtual points strategy is used to simplify the model. The mathematical model is established based on VRPMTW (Favaretto et al., 2007). Let G = ( , ) be a graph. Let V = 0 ∪ N be the set of distribution center and work stations. Let = (V × ) be the set of arcs: there are travel times associated with each arc. Let D = , , … be the demands set of stations. If the demands are expressed in the minimum distribution unit, then D = , , , … , , , … , … , , … . Let N + be the set of stations that require multiple visits. Station i of this set has demand time windows, [ , ], 1 ≤ v ≤ . The time windows of a station do not overlap with each other: < < , v ∈ {1,2, … , − 1}. This time interval is not less than the ratio of the planning period to the number of packages + 1 ⁄ . Each station with N + is provided virtual points that are associated with real stations and real time windows [ , ]. They are in the same location. All the stations of M must be served once.
The objective function (5) is to minimize the sum of the fixed costs of the vehicles, the waiting time costs and the transportation time costs. To convert the cost of using a vehicle to the same units used in the other parts, the fixed cost is expressed in time units.
Constraints (6)-(7) state that each virtual station is assigned to exactly only one vehicle. Constraint (8) states that the number of arcs leaving station i is equal to the number of arcs entering it. Constraints (9)-(10) state that each vehicle used for distribution starts and ends at the depot. Constraint (11) states that the maximum amount of loaded material is upper bounded by the capacity of vehicle type traversing arc(i, j). Constraints (12)-(13) mean that the arrival time at station is equal to the arrival time at station plus the waiting time and service time plus the travel time, only if this arc is assigned to vehicle . Constraints (14)-(15) mean that the arrival time at station plus the waiting time must meet the time window of station . Constraint (16) means that the number of vehicles must not exceed the prescribed number. Constraints (17)-(19) state the feasibility intervals for the decision variables.
The mathematical model describes the problem clearly and intuitively by using the virtual points strategy. There are many features implicit in the model that are not intuitively expressed. Constraints (20)-(21) state that the real stations may need multiple visits, and the vehicle entry and exit times are equal to the number of time windows. Constraint (22) means that the demand of the stations visited by more than one vehicle is satisfied. Constraints (23) means that the material delivered by vehicle must be less than the stock limit of the line. The linear solver can solve this problem small-scale in a short computational time. However, for the practical production environment, the size of the proposed logistics planning problem is relatively large. And the problem may need to be solved and updated several times in each production horizon (e.g. a day). It is difficult to obtain the solutions in a reasonable computational time for the linear solver. Thus, a meta-heuristic algorithm is designed to address such an operational problem with fast speed for the real application scene. Dantzig and Ramser (1959) proved that the VRP is an NP-hard problem. Therefore, SHVRPMTW is also an NP-hard problem but more complex. Improved meta-heuristic algorithms are effective methods to address logistics planning problems. Among these algorithms, ACO is especially outstanding: ACO models a group of ants that cooperate through information exchange by iterating pheromones on the edge of a graph. Based on the OVPR algorithm of Li et al. (2009), HACO-TS is improved to solve the problem effectively.

Hybrid Ant Colony Optimization Algorithm Based on Tabu Search
The transfer rules are improved according to the characteristics of the problem. The solutions that meet the targets are quickly identified under the surrounding search strategy. The short-term memory capability and defiance criteria of the TS are incorporated into the algorithm. HACO-TS makes good use of TS's strong local search ability and ACO's parallel global search ability to improve the convergence performance. HACO-TS algorithm process is shown in Fig. 3. Jiang et al. (2014) proposed that the service order of customers and vehicles affects the quality of the initial solution. They classified customers and vehicles into four categories. To minimize the time cost, the greatest distance rule can be adopted.

Sort stations
Referring to the vehicle classification rules, the greatest capacity rule is adopted, meaning that larger vehicles are preferred.

Algorithm parameter initialization
First station is selected  Fig. 3. HACO-TS algorithm flow chart

4.1.1Surrounding search strategy
The feasible solution is preliminarily constructed by combining the minimum total cost of the target and vehicle sequencing rules. Let ℶ be the feasible points.
When solving large-scale problems, it is difficult to search for the optimal solution because the algorithm faces a large search space. The surrounding search strategy is improved to narrow the search scope. All points in the set of the next feasible points are sorted according to the Euclidean distance, and points in a certain range are selected as the peripheral point set ( ), which is the new candidate points set. The neighborhood operators are modified according to this strategy. The concept of geographic proximity (Br Ysy et al., 2009) is defined, which represents the average X and Y coordinates of the stations on each route.

4.1.2First station selection strategy
To select a suitable node from many feasible nodes, the first station selection strategy is used to improve the convergence of the algorithm. When the iteration progress reaches a certain condition, the first node is the node with the highest transition probability in the candidate node set. The strategy is defined as follows: where is the first station of each route.
(t) is the probability from the distribution center to using vehicle K. is the number of current iterations.
is the maximum number of iterations. ε(0 < ε < 1) indicates when to use the strategy.

Solution construction
After all the ants have completed the itinerant routes, the pheromone updating rule is used to update the pheromone level, thus increasing the pheromone on the edge belonging to the current optimal itinerant routes. HACO-TS ends when one of the following conditions is met: (1) a fixed number of solutions are generated, (2) a fixed number of iterations have not executed the objective function.

Probabilistic transition rule
Each ant selects stations according to the probabilistic transition rule to generate a complete itinerant route. To enhance the accuracy and reality of ants' route selection, the transition rule is improved on the basis of the adaptive pseudorandom ratio selection rule proposed by Li et al. (2009).
is the pheromone concentration. , and indicate the relative importance of pheromones and heuristics.
(0 ≤ ≤ 1) is a constant to control the transition rule. q is a random variable subject to a standard uniform distribution.
When q ≤ , the station corresponding to the maximum probability is selected as the next visit station. Otherwise, station is selected. is obtained by roulette selection, where the probability is (t). is the expected heuristic function, representing the expected degree of the ant's transfer from station i to station j.
takes into account not only the time but also the urgency of station j. The urgency is given by the time interval between the current moment and the current time window of station j.
where represents the current time starting from station i and is a scale factor with uniform mass. ( ) = 1 ( , ) ⁄ , ∈ ( ) represents a function related to wait(i, j) of the vehicle at station j.

The pheromone update strategy
In each iteration of the algorithm, the pheromone matrix is updated after all the ants have constructed the solution. The pheromone update strategy: The increment is calculated by mixing the best solution found in the current iteration and the global optimal solution , which is improved based on Li et al.(2009). The pheromone delta ∆ emitted by ants: where, is the minimum objective value corresponding to , and is the minimum objective value corresponding to . and are weight parameters, representing the relative importance of the solution in pheromones increment, > 0, > 0, + = 1.When = 0, = 1 indicating only is used to update the pheromone, the convergence speed of the algorithm is the fastest, but it is easy to fall into local optima. Only if ( , ) ∈ , ( , ( , ) = 1; else, ( , ( , ) = 0. ρ(0 < ρ ≤1) denotes the volatility coefficient of pheromones on the routes, and 1 − ρ denotes the persistence of pheromones. In the process of updating pheromones, the pheromone on each edge ( , ) evaporates at a fixed evaporation rate. Different combinations of and determine the different updating trajectories of the arc pheromones. As the iteration process progresses, the importance of decreases and the importance of increases; therefore, the relative weight coefficients are adjusted dynamically.
Based on the MMAS pheromone, the pheromone is controlled between [ , ] to prevent premature convergence of the algorithm (Stützle & Hoos, 2000). The pheromone needs to be judged after each iteration. = * , = .

Local search
Numerous studies have shown that the ACO integrated local search can lead to better solutions (Li et al., 2009, Liu and. TS is an efficient search algorithm to solve the VRP (Taillard et al., 1997;Xia & Fu, 2018a, b;Gendreau et al., 1999); TS with multiple starts is used as the local search here. A good balance between the quality and the computation time must be considered. If all the solutions found by ants are searched locally, the process will require considerable computation time to obtain better results. The number of ants in the local search gradually increases with every iteration base on the strategy of Referring to the split operation of Jiang et al. (2014) to optimize the vehicle types used. The chosen route is randomly divided into two parts. The algorithm takes the FA standard if the split reduces the total cost. During the operation, the work order of the route remains the same. To avoid circuitous search, the tabu list must be set to prohibit the solution that was just selected. The tabu objects are vertex pairs in operations. To increase the randomness, the tabu length is set as a random integer in the interval [ , ]. The tabu lists adopt the first-in-first-out strategy. To avoid excessive tabu, a tabu list reinitialization strategy is designed. The evaluation function is implemented by a hierarchical calculation: the fixed cost of the vehicle is given priority, followed by transportation time and waiting time. The key point of local search is to further optimize the solution in less computation time. The standard tabu search stores the preferred solution, while this article's local search stores attributes.

Postprocessing
To further optimize the optimal solution found by the algorithm without significantly affecting the computation time, the algorithm integrates a postprocessing process. VNS has been successfully applied to VRP variants as an intensive search strategy. The process is based on an improved VNS reinforcement strategy. When the main loop terminates, postprocessing is implemented to improve the optimal solution. The best solution is used as the initial solution of VNS. The local optimal solution of one kind of neighborhood operator is not necessarily the optimal solution of another one. The local optimal solution under all neighborhood operators should be the global optimal solution. VNS includes two main stages: shake and local search. The shake mechanism is used as a diversification strategy to further explore the search space. The goal of postprocessing is to reduce the total distance and time rather than the number of vehicles. The key point of postprocessing is to further optimize the solution in less computation time.

The time relaxation rule
The time relaxation rule is designed based on the problem characteristic that the time window cannot be violated. This rule is used in decoding to update the minimum wait time and adjust the arrival time to the stations in the routes in a backwards manner. The service time of the stations is within the time window and conforms to the constraints. The delay of the arrival time ensures the latest leave time of vehicles and does not violate the time window.

Algorithm testing
This section describes the analysis of the results to evaluate the performance of the algorithm. The experiments are conducted in Java, on a computer with Intel Core (TM) i5-7500 CPU @3.40GHz and 8 GB RAM.

Testing problems
This proposed problem is a deformation of VRPTW in two dimensions of vehicle types and time windows. There are currently no standard test cases for this problem. The test cases are generated based on Gendreau et al. (1999) and Jiang et al. (2014). Due to the high fixed cost of the vehicles in the problem, C4, C6, and C14-C20 are selected and deformed. The locations of stations, total demand and service time are all derived from the cases. The test cases are shown in Appendix1. The second time split has more time windows than the first time.

Parameter settings
Because the parameters are very important to the performance of the algorithm, this section designs the Taguchi experiment to choose the best combination of parameters. The following parameters are controlled: , , q , ρ . Each parameter combination is tested 10 times. For all experiments, the relative percentage deviation (RPD) was used as the response variable.
where is the best value of the same instance in the parameter tests and is the value of the function obtained in the experiment. It can be seen that the lower RPD, the better the performance of parameter combination. As shown in Table.1, there is a range of parameters and (4 )parameter combinations. The range of parameters are: ∈ {2,3,4,5}, ∈ {1,2,3,4}, q ∈ {0.6,0.8,0.9,0.95}, ρ ∈ {0.05,0.1,0.2,0.35}. The effect diagram of the result analysis is shown in Fig.6. The signal-to-noise ratio (SNR) is the ratio of the objective function value to the variance of the objective function. By comparing the slopes of the lines, the relative magnitude of each factor's influence can be assessed. The number of shared solutions has the greatest impact on the optimal solution. The optimal parameter combination of the current algorithm: = 3, = 1, q = 0.8, ρ = 0.2. Fig. 6. The average main effects plot and SNR main effects plot

Comparisons with CPLEX
The mathematical model is using IBM ILOG CPLEX12.5 to evaluate the accuracy. The results calculated by CPLEX and the improved algorithm are compared, as shown in Table 2 although the result of the linear solution is the best. For normal-scale problems, only the solving speed of the meta-heuristics algorithm is likely to be acceptable. The advantage of the computing cost for meta-heuristics is more obvious than CPLEX. The hybrid strategy can be compared with only the contemporary optimal solution and the global optimal solution. Note that none of these approaches use postprocessing. In HACO-TS-sb, = 0 means that only the global optimal solution is considered to update the pheromone, and in HACO-TS-ib means = 1. HACO-TS-wpp uses the hybrid strategy to update pheromones. The averaged results are recorded in Table 3. The last line gives the ratio of the algorithm to the other algorithms.

Table 3
Comparison of the pheromone update strategy Considering the balance between the solution improvement and the computation time, adding the post optimization is of great significance for the further improvement of the solution. Fig.7 shows Tukey's HSD interval plots at 95% confidence level for the algorithms. Because the fixed cost is included in the result statistics, the APRD difference is small. It is clear that the ARPD of HACO-TS is lower than other variants. This indicates that HACO-TS is the closest to the optimal value, and there is a significant statistical difference. The ARPD of HACO-TS is 0.0013, the other variants are respectively 0.0058, 0.0066, 0.0050. This finding proves statistically that the postprocessing and neighborhood operations could significantly improve the performance of HACO-TS.

Effectiveness of HACO-TS
All the selected algorithms are used to solve the heterogeneous VRP with one constraint loosened: the same station can be visited multiple times. The algorithms are as follows: SA is a deterministic annealing meta-heuristic proposed by Bräysy et al.(2008); VNTS is a reactive variable neighborhood tabu search proposed by Paraskevopoulos et al.(2008); TTS is the two-phase tabu-search algorithm proposed by Jiang et al.(2014).
To ensure fairness, these algorithms use the same goals and the number of iterations for comparison. Each case is executed 10 times, and the statistical optimums are shown in Table 4. Comparison with similar algorithms shows that HACO-TS is superior and effective. The comparison of 1 and 2 for each case shows that the more containers in the same instance, the higher delivery costs because longer routes or more waiting time are required to meet the demand. Even under the limitation of time windows, some instances such as C19-1 and C19-2, need to increase the number of vehicles to meet the requirements.
Specifically, the following conclusions can be drawn: (1) Based on the goal orientation, no algorithm outperforms all other algorithms in all cases. However, from the perspective of the optimal solution overall, the solution quality of HACO-TS is high.
(2) For most cases, for example C17-C20, the optimal solution found by the HACO-TS is the best. In general, the optimal solution requires fewer vehicles.
At the beginning, the initial setting of the pheromone and the first station selection strategy can improve the global search capability of the algorithm. The probabilistic transition rule and the pheromone update strategy can effectively improve the effectiveness of the algorithm. The pheromone update strategy greatly improves the diversity of the algorithm search process and promotes the exploration of the solution space. The quality of the solution is improved by local search, which illustrates the good performance of HACO-TS.   Fig. 8 shows that HACO-TS produces the optimums for almost all cases and is closest to the optimal APRD value. It is further proved that HACO-TS is better than other algorithms. Different algorithms use different programming languages and on different machines, which greatly affects the calculation time. It is difficult to fairly compare the computation times of different algorithms. Differences in programming ability can affect the calculation time. Moreover, the programming language, algorithmic mechanisms, and data structures all impact the computation time. HACO-TS is a meta-heuristics algorithm with high accuracy that combines the advantages of ACO and TS. The effectiveness of HACO-TS is verified by comparing the results. The better solution of ACO is set as the initial solution of postprocessing, which increases the possibility to obtain a better solution. However, the computation time is relatively burdensome, and the required accuracy of the question is important to consider. Compared with the computation time, the performance of the algorithm is more consistent with the evaluation index.

Case study
This section uses data of an enterprise for a case study. The rationality of the plan is the foundation of production. The total demands of stations can be calculated according to the production plan. There are 1 center and 25 stations that make supplementary material requests in the workshop. There are 3 types of vehicles shown in Table 5. Based on inventory management, the upper limit for each line-side is obtained. According to the material properties and production plan, the consumption time of each unit of material and multiple call times are calculated. The distribution plan is configured according to the expected time windows. Vehicles in the mixed-flow manufacturing workshop need to complete multiple distribution tasks. In contrast to the general VRP, how to complete the distribution route assignment with the least cost must be considered.  B 0-14-16-18(1)-8-12(2)-6(2)-0 2 A 0-13(1)-21(1)-15-22-13(2)-0 3 A 0-5-2(1)-7-10-20-1(2)-0 4 A 0-12(1)-1(1)-9-3-21(2)-2(2)-0 5 A 0-6(1)-11-19-17-18(2)-0 6 A 0-23-4-24-25-0 The same parameter configuration is designed by Taguchi experiment. The results are shown in Table.6, Table.7. The number of vehicles is gradually reduced as the vehicle load increases. Considering the high cost of large-capacity vehicles, hybrid types are adopted for delivery to reduce the fixed cost. A satisfactory solution can be obtained if the number of distribution vehicles remains constant. In the actual distribution, considering the speed, low energy consumption of small vehicles and increasingly high demands on the timeliness of distribution, the multivehicle combination distribution scheme will improve the service quality and reduce the distribution cost to a large extent. This method produces reasonable planning for the distribution and achieves efficient logistics distribution. The research has greatly reduced delivery costs by delivering materials in the least amount of time to meet production needs. Moreover, the research can be used in actual production environments.

Conclusion
Smooth material flow is an important condition for enterprises to realize just-in-time production. SHVRPMTW is a realistic VRP with wide application value in manufacturing systems. Combined with the material requirement characteristics, the article studies SHVRPMTW for the first time and establish the mathematical model. Modeling with the virtual point strategy reduces the complexity of the problem. HACO-TS is designed to solve large-scale problems. The multistart strategy enhances the diversity and increases the likelihood of finding better solutions. Local search and postprocessing are performed by combining several neighborhood operations. Based on the problem characteristic that the time window cannot be violated, the time relaxation rule is designed to update the minimum wait time. Based on test cases, the effectiveness of HACO-TS is verified by comparing the results with other algorithms.
In future research, various uncertain factors, such as the uncertain variety and quantity of demands, random arrival time, and even the abnormal vehicle transportation, can be added to the vehicle route planning problem with split deliveries for station demands.

Statements and Declarations
We declare that we have no competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.