Medicine logistics vehicle path optimization method based on two-stage heuristic algorithm

Abstract

The invention relates to a medicine logistics vehicle path optimization method based on a two-stage heuristic algorithm, and belongs to the field of vehicle path optimization. The method comprises the following steps of: collecting longitude and latitude information of a warehouse and a client point on a specific route of a certain medical company by using a map API (Application Program Interface) so as to calculate a distance matrix between the client point and the warehouse; and establishing a medicine logistics vehicle path optimization model with time window limitation, wherein the optimization target is to minimize the transportation cost and the time window penalty cost of the vehicle. A two-stage heuristic algorithm is designed, in the first stage, on the basis of considering client point time windows and geographic positions, client space-time distances are calculated, and clustering partitioning is conducted on the client space-time distances through an improved K-means clustering algorithm; in the second stage, a genetic algorithm improved by a self-adaptive large-field search algorithm is combined to optimize the transportation route of each region; the objective function of the medical logistics vehicle path optimization model is solved. According to the invention, an efficient logistics driving route is generated, and the cost is reduced.

Description

Technical Field

[0001] The invention belongs to the field of vehicle path optimization and relates to a pharmaceutical logistics vehicle path optimization method based on a two-stage heuristic algorithm. Background Art

[0002] At present, optimizing the distribution route has become a key link in improving logistics efficiency and reducing transportation costs. Effective route selection can not only improve distribution efficiency, but also ensure the timely supply of drugs, improve customer satisfaction, and enhance the market competitiveness of pharmaceutical companies. Therefore, the development of advanced vehicle route optimization technology will provide important support for the sustainable development of pharmaceutical logistics and help the industry achieve healthy development in the fierce competition.

[0003] In the field of vehicle routing optimization, heuristic algorithms are often used to solve complex vehicle routing optimization problems, especially for NP-hard problems. Since the time cost of classical algorithms in solving such problems is too high and often difficult to meet practical needs, it is crucial to adopt more effective solution methods. Heuristic algorithms provide a strategy to approximate the optimal solution within an acceptable computing time. Although it cannot ensure that the optimal solution is found, in many practical application scenarios, the approximate solution obtained is often sufficient to meet the needs.

[0004] Heuristic algorithms usually start from a feasible solution and repeatedly search for a better solution to approach the optimal solution within a reasonable computing time. However, for large-scale NP-hard problems, if the traditional heuristic algorithm is not improved, the quality of the solution and the efficiency of the solution may not be ideal. Patent CN110852469A discloses a method for optimizing the design of pharmaceutical logistics distribution paths based on a genetic algorithm, which includes constructing a parallel genetic algorithm framework based on a MapReduce model for the pharmaceutical logistics distribution path problem. However, this method still has the problem that the feasible solution found at the beginning may be poor, resulting in a poor approximate optimal solution obtained in the end. Therefore, the present invention proposes a pharmaceutical logistics distribution path optimization method based on a two-stage heuristic algorithm. In the first stage, an improved k-means clustering algorithm is used to effectively cluster customer points, thereby reducing the scale of the problem; in the second stage, a genetic algorithm improved by combining an adaptive large field search algorithm is used for further optimization. This method effectively solves the problem of poor quality of the initial feasible solution, improves the quality of the approximate optimal solution finally obtained, and thus improves the overall efficiency of pharmaceutical logistics distribution. Summary of the invention

[0005] In view of this, the purpose of the present invention is to provide a pharmaceutical logistics vehicle path optimization method based on a two-stage heuristic algorithm, which effectively classifies customers with similar spatial and temporal distances into one category, and then further optimizes each area in combination with a genetic algorithm improved by an adaptive large-area search algorithm. This vehicle path optimization method based on the idea of ​​"clustering first, then optimizing" can effectively improve the results of pharmaceutical distribution path optimization. The present invention can effectively reduce distribution costs, improve distribution efficiency, and enhance customer satisfaction.

[0006] In order to achieve the above object, the present invention provides the following technical solutions:

[0007] A pharmaceutical logistics vehicle routing optimization method based on a two-stage heuristic algorithm "clustering first, then optimizing" for all customers can effectively search for the distribution plan with the lowest total cost in the entire area.

[0008] The first stage of the improved K-means clustering algorithm includes the following steps:

[0009] S1: Use the map API (Application Programming Interface) to collect the longitude and latitude information of a certain route of a pharmaceutical company, including warehouses and all customer points, and use the map API to convert the collected longitude and latitude information of warehouses and all customer points into metric coordinates.

[0010] S2: Randomly select K points as the initial cluster centers, where represents the total demand of customers, Q max Indicates the maximum load weight of the vehicle.

[0011] S3: Calculate the spatiotemporal distances between the remaining customer locations and the K cluster centers.

[0012] Since the dimensions of time distance and space distance are different, time and space must be normalized before calculating time and space distance. After normalization, the time and space distance between two customer points can be calculated by weighted average, where α is the weight update coefficient of time distance and space distance.

[0013]

[0014] The spatial distance between customer points can be described by the Euclidean distance. Assume that the coordinates of customer point i and customer point j in the two-dimensional plane are (x i ,y i ) and (x j ,y j )Then the Euclidean distance between customer point i and customer point j is:

[0015]

[0016] The time distance between customer point i and customer point j is described by the satisfaction degree of the service time window. A soft time window strategy is adopted to introduce a certain penalty coefficient for the time distance beyond the time window range to characterize it. The time windows of customer point i and customer point j are [e i ,l i ] and [e j ,l j ], let e i <e j , if a i ∈[e i ,l i ], then the time a that the vehicle arrives at customer point j is j ∈[e i +s i +t ij ,l i +s i +t ij ], remember e j ′=e i +s i +t ij ,l j ′=l i +s i +t ij , the penalty coefficient for lateness is defined to be greater than the penalty coefficient for waiting (τ2>τ1), then the time distance between customer points i and j is:

[0017]

[0018] S4: Based on the vehicle load constraint, the remaining customers are assigned to the partition corresponding to the nearest centroid according to the time and space distance, and S5 is executed for customer points that exceed the vehicle load constraint;

[0019] S5: Allocate customers who do not meet the current load constraint to other partitions with the smallest spatiotemporal distance from them.

[0020] S6: Check whether there are customers who have not been delivered. If so, jump to S4, otherwise recalculate the centroids of K partitions and repeat the above steps until the clustering results converge.

[0021] The K-means clustering algorithm is used to cluster the customer set based on the time-space distance, laying the foundation for the subsequent delivery route optimization.

[0022] The second stage combines the improved genetic algorithm with the adaptive large domain search algorithm and includes the following steps:

[0023] S1: Use the map API to calculate the distance matrix between the warehouse and all customer points’ longitude and latitude information collected by the algorithm in the previous stage.

[0024] Using the map API and the longitude and latitude information of the queried warehouse and customer points, calculate the round-trip distance between every two points.

[0025] The queried round-trip distance is stored in a (warehouse + all customer points)^2 matrix.

[0026] Export the matrix as a CSV file, which is the distance matrix containing warehouses and all customer points on a certain route of a pharmaceutical company.

[0027] S2: Construct a pharmaceutical logistics vehicle path optimization model, with the optimization goal of minimizing the transportation cost and time window penalty cost. The vehicle's driving route is the total distance from the warehouse, traversing all customer points, and finally returning to the warehouse; the vehicle has a capacity limit and no maximum distance limit, and set constraints on the objective function of the pharmaceutical logistics vehicle path optimization model with time window restrictions.

[0028] The overall optimization objective function of the pharmaceutical logistics vehicle routing optimization model with time window constraints described in S1 is:

[0029]

[0030] Where N = {0, 1, 2, ... n} is the set of distribution centers and customer points, N′ = {1, 2, 3, ... n} is the set of customer points, K = {1, 2, ... k} is the set of available vehicles in the warehouse, i, j represents a single warehouse or customer point, i, j ∈ N, i ≠ j; represents the unit distance cost of vehicle k traveling from customer point i to customer point j, d ij represents the actual distance between node i and node j, α1 represents the penalty cost of early vehicle arrival, and α2 represents the penalty cost of late vehicle arrival; e i represents the left time window of customer i, t i represents the time when the vehicle arrives at customer i, l i represents the right time window of customer i; is a decision variable 0-1 variable, which takes 1 when vehicle k is traveling on path (i, j), otherwise it takes 0

[0031] The objective function setting constraints of the pharmaceutical logistics vehicle routing optimization model with time window constraints described in S2 are as follows:

[0032]

[0034] τ0=e0 (6)

[0035]

[0036] Formula (1) indicates that the distribution center has enough vehicles to complete the distribution task; Formula (2) indicates that every customer is served; Formula (3) indicates that each vehicle starts from the distribution center and finally returns to the distribution center; Formula (4) indicates that the total demand of customers served by each vehicle must not exceed the maximum load of the vehicle; Formula (5) indicates the demand of customer point i is met; Formula (6) indicates that the departure time of the vehicle is the earliest opening time of the distribution center; Formula (7) indicates that the time when the vehicle returns to the distribution center must be earlier than the closing time of the distribution center; Formula (8) indicates that the time when the vehicle reaches the next node is equal to the sum of the time when it arrives at the previous node, the service time of the previous node and the travel time between the two nodes; Formula (9) is the value constraint of the decision variable.

[0037] S3: Design a genetic algorithm combined with an adaptive large neighborhood search algorithm. In the first stage, the greedy algorithm is used to generate the initial solution. In the second stage, the genetic algorithm is used to optimize the generated initial solution. In the third stage, the adaptive large neighborhood search algorithm is used to perform deep optimization. The solution optimized by the genetic algorithm is used to avoid falling into the local optimum of the genetic algorithm to obtain the optimal distribution plan. The genetic algorithm (Genetic Algorithm AdaptiveLarge Neighborhood Search, GA-ALNS) improved by the adaptive large neighborhood search algorithm is used to solve the objective function of the pharmaceutical logistics vehicle path optimization with time window constraints:

[0038] The specific steps of the GA-ALNS algorithm mentioned in S3 are as follows:

[0039] Step 1: Use the greedy algorithm to generate the initial solution: select a customer point i that is closest to the distribution center in the customer set, add it to the individual, then search among the customers who have not yet been added to the individual, find the customer j that is closest to the current customer and add it to the individual, and repeat this operation until all customers are added to the individual.

[0040] Since the transport vehicle has a weight limit, the algorithm should follow the following principles when generating the initial population: when the total demand of customers is less than the maximum load of the vehicle when the number of customers is accumulated to n+1, the result is the opposite. In this case, the number of codes of the distribution center is inserted after the nth customer, indicating that the coding of a chromosome is completed. For example, if there are 5 customers in total, after the first and second customers are inserted, the third customer is inserted, which violates the weight constraint. In this case, the chromosome code is 1-2-6, and this is repeated until all individuals are included to form a number that can meet the population size.

[0041] Step 2: Calculate the fitness value of the population. The higher the fitness value of an individual in the population, the better the individual's performance. Since the objective function set in this paper is to find the minimum value, the fitness function is taken as the inverse of the objective function.

[0042] Step 3: Perform selection, crossover, and mutation operations on all individuals in the population to obtain a new population P1. The selection operation here is a binary tournament selection operation, the crossover operation is an OX crossover method, and the mutation operation is a two-point mutation.

[0043] Step 4: Randomly select R chromosomes from the new population P1 as another population P1.1, and the remaining chromosomes are recorded as P1.2.

[0044] Step 5: Population P1.1 is recorded as Li and enters the Adaptive Large Neighborhood Search (ALNS) operation. First, the ALNS parameters are initialized, and then all chromosomes in Li are destroyed and reorganized through the destruction operator and repair operator in ALNS.

[0045] Remove the operator set h r as follows:

[0046] (1) Similarity removal operator: remove the node with the smallest correlation with the selected node, that is, randomly select a customer i to remove, calculate the correlation D between the remaining customers and i, and select the customer with the largest correlation to remove, that is,

[0047]

[0048] Where V ij Is i and j on the same road, that is, whether they are served by the same vehicle. If i and j are on the same road, it is 1, otherwise it is 0; d in Represents the distance from customer i to any customer.

[0049] (2) Random loop removal operator: selects a path with the largest total cost and then randomly selects m customers in the path to remove.

[0050] (3) Worst distance removal operator: Calculate the sum of the distances between each customer and the customers before and after it on the vehicle path, and select the top m customers with the largest sum of distances.

[0051] (4) Worst time removal operator: select the customer with the largest time deviation and remove it from the current path. The time deviation of customer j WT(j) = |y j -a j |, where y j is the time when the delivery vehicle arrives at customer j, a j is the earliest acceptable service time for customer j.

[0052] Insert operator set h i as follows:

[0053] (1) Optimal greedy insertion operator: Find feasible insertion locations for the destroyed customer points one by one, that is, satisfy the time window and load constraints, and select the location with the smallest path cost increment among all feasible insertion locations for insertion.

[0054] (2) Random insertion operator: randomly insert any position of the solution after destruction under the premise of satisfying the time window and load constraints.

[0055] (3) Optimal distance insertion: Find feasible insertion locations for the destroyed customer points one by one, that is, satisfy the time window and load constraints, and select the location with the smallest distance increment from all feasible insertion locations for insertion.

[0056] (4) Best time insertion heuristic operator: For each subpath in the solution to be reconstructed, under the constraints of time window and capacity, first find all feasible insertion positions for customer j and calculate the time deviation after insertion. The time deviation of customer j, BT(j) = |y j -a j |, where: y j is the time when the delivery vehicle arrives at customer j, a j is the earliest acceptable service time for customer j, and then customer j is inserted at the position with the smallest time deviation.

[0057] Step 5: When performing repair and destruction, each time a roulette wheel is used to select one operator from the four destruction operators and four repair operators for destruction and repair. Therefore, an adaptive mechanism is designed to evaluate the effects of operators by scoring and change the probability of selecting them. The specific operation is as follows:

[0058] The weight of the current operator participating in the roulette selection is calculated based on the average score of each round. Initially, the weight of each operator is set to 1, and the initial score is set to 5. If the solution found by the operator is better than the current global optimal solution, the operator will add α1 points; if the solution found by the operator is better than the current solution, the operator will add α2 points; if the solution found by the operator is worse than the current solution, the operator will add α3 points. In addition, a weight update coefficient ρ is set to avoid falling into the local optimum due to too fast convergence speed. The weight of the operator is updated according to the following formula:

[0059]

[0060] Where: a is the weight, s a is the cumulative score, μ a is the cumulative number of times used, and ρ is the weight update coefficient in the adaptive mechanism.

[0061] Step 6: Use the following rules to control the update of the solution after ALNS. Let s1 represent the current solution, s2 represent the new solution, f(s1) and f(s2) be the objective function values ​​of s1 and s2 respectively. If f(s2)-f(s1)<0, replace the current solution with the new solution, otherwise reject the current solution. ALNS ends when all chromosomes in Li have been ALNSed.

[0062] Step 7: Reorganize P1.1 and Li after ALNS operation into a population P1. Because in the selection operation, one duplicate individual is retained and the rest are deleted. Therefore, after crossover, mutation, and ALNS operations, the number of new offspring is less than the original population. Recombination operation is required to combine P1 with the previous generation population to keep the population unchanged.

[0063] Step 8: Determine whether the maximum number of iterations has been reached. If not, return to step 2 until the maximum number of iterations has been reached.

[0064] The beneficial effects of the present invention are:

[0065] (1) The present invention uses a two-stage heuristic algorithm of "clustering first, then optimizing" to effectively classify customers with similar spatial and temporal distances into one category, thereby reducing the problem scale and improving the solution efficiency.

[0066] (2) The present invention adopts an improved K-means clustering algorithm and utilizes time-space distance to cluster customer points. This can more accurately reflect the actual needs of customers than the traditional clustering algorithm based on spatial distance, thereby improving the quality of clustering results.

[0067] (3) The present invention combines the genetic algorithm improved by the adaptive large-domain search algorithm, effectively avoiding the genetic algorithm from falling into the local optimal solution, improving the quality of the final approximate optimal solution, and thus improving the overall efficiency of pharmaceutical logistics distribution.

[0068] (4) The present invention can effectively reduce distribution costs, improve distribution efficiency, and enhance customer satisfaction, providing important support for the sustainable development of pharmaceutical logistics.

[0069] Other advantages, objectives and features of the present invention will be described in the following description to some extent, and to some extent, will be obvious to those skilled in the art based on the following examination and study, or can be taught from the practice of the present invention. The objectives and other advantages of the present invention can be realized and obtained through the following description. BRIEF DESCRIPTION OF THE DRAWINGS

[0070] In order to make the purpose, technical solutions and advantages of the present invention more clear, the present invention will be described in detail below in conjunction with the accompanying drawings, wherein:

[0071] Figure 1 Flowchart of the K-means clustering algorithm based on spatiotemporal clustering.

[0072] Figure 2 Flowchart of the improved genetic algorithm combined with the adaptive large domain search algorithm.

[0073] Figure 3 Flowchart for generating initial solutions for the greedy algorithm. DETAILED DESCRIPTION

[0074] The following describes the embodiments of the present invention by specific examples, and those skilled in the art can easily understand other advantages and effects of the present invention from the contents disclosed in this specification. The present invention can also be implemented or applied through other different specific embodiments, and the details in this specification can also be modified or changed in various ways based on different viewpoints and applications without departing from the spirit of the present invention. It should be noted that the illustrations provided in the following embodiments only illustrate the basic concept of the present invention in a schematic manner, and the following embodiments and features in the embodiments can be combined with each other without conflict.

[0075] Among them, the drawings are only used for illustrative explanations, and they only represent schematic diagrams rather than actual pictures, and should not be understood as limitations on the present invention. In order to better illustrate the embodiments of the present invention, some parts of the drawings may be omitted, enlarged or reduced, and do not represent the size of actual products. For those skilled in the art, it is understandable that some well-known structures and their descriptions in the drawings may be omitted.

[0076] The same or similar numbers in the drawings of the embodiments of the present invention correspond to the same or similar parts; in the description of the present invention, it should be understood that if the terms "upper", "lower", "left", "right", "front", "rear", etc. indicate the orientation or position relationship, they are based on the orientation or position relationship shown in the drawings, which is only for the convenience of describing the present invention and simplifying the description, rather than indicating or implying that the device or element referred to must have a specific orientation, be constructed and operate in a specific orientation. Therefore, the terms describing the position relationship in the drawings are only used for illustrative purposes and cannot be understood as limiting the present invention. For ordinary technicians in this field, the specific meanings of the above terms can be understood according to specific circumstances.

[0077] The technical solution adopted by the present invention to solve the technical problem is: a method for optimizing the path of pharmaceutical logistics vehicles based on a two-stage heuristic algorithm, comprising the following steps:

[0078] The first stage algorithm:

[0079] S1: Use the map API to collect the longitude and latitude information of a certain route of a pharmaceutical company, including warehouses and all customer points, and use the map API to convert the collected longitude and latitude information of warehouses and all customer points into metric coordinates.

[0080] S2: Randomly select K points as the initial cluster centers, where represents the total demand of customers, Q max Indicates the maximum load weight of the vehicle.

[0081] S3: Calculate the spatiotemporal distances between the remaining customer locations and the K cluster centers.

[0082] S31 Since the dimensions of time distance and space distance are different, time and space must be normalized before calculating time and space distance. After normalization, the time and space distance between two customer points can be calculated by weighted average.

[0083]

[0084] The spatial distance between S32 customer points can be described by the Euclidean distance. Assume that the coordinates of customer point i and customer point j in the two-dimensional plane are (x i ,y i ) and (x j ,y j )Then the Euclidean distance between customer point i and customer point j is:

[0085]

[0086] The time distance between the third customer point i and customer point j in S3 is described by the satisfaction degree of the service time window. A soft time window strategy is adopted to introduce a certain penalty coefficient for the time distance beyond the time window range to characterize it. The time windows of customer points i and j are [e i ,l i ] and [e j ,l j ], let e i <e j , if a i ∈[e i ,l i ], then the time a that the vehicle arrives at customer point j is j ∈[e i +s i +t ij ,l i +s i +t ij ], remember e j ′=e i +s i +tij ,l j ′=l i +s i +t ij , the penalty coefficient for lateness is defined to be greater than the penalty coefficient for waiting (τ2>τ1), then the time distance between customer points i and j is:

[0087]

[0088] S4: Based on the vehicle load constraint, the remaining customers are assigned to the partition corresponding to the nearest centroid according to the time and space distance, and S5 is executed for customer points that exceed the vehicle load constraint;

[0089] S5: Allocate customers who do not meet the current load constraint to other partitions with the smallest spatiotemporal distance from them.

[0090] S6: Check whether there are customers who have not been delivered. If so, jump to S4, otherwise recalculate the centroids of K partitions and repeat the above steps until the clustering results converge.

[0091] Second stage algorithm:

[0092] S7: Use the map API to calculate the distance matrix between the warehouse and all customer points’ longitude and latitude information collected by the algorithm in the previous stage.

[0093] S71 uses the map API and the longitude and latitude information of the queried warehouse and customer points to calculate the round-trip distance between every two points.

[0094] S72 stores the queried round-trip distance in a (warehouse + all customer points)^2 matrix.

[0095] S73 exports the matrix as a CSV file, which is a distance matrix including warehouses and all customer points on a certain route of a pharmaceutical company.

[0096] S8: Construct a pharmaceutical logistics vehicle path optimization model, with the optimization objective of minimizing the transportation cost and time window penalty cost. The vehicle's driving route is the total distance from the warehouse, traversing all customer points, and finally returning to the warehouse; the vehicle has a capacity limit and no maximum distance limit, and set constraints on the objective function of the pharmaceutical logistics vehicle path optimization model with time window restrictions.

[0097] S81: The overall optimization objective function of the pharmaceutical logistics vehicle routing optimization model with time window constraints described in S8 is:

[0098]

[0099] Where N = {0, 1, 2, ... n} is the set of distribution centers and customer points, N′ = {1, 2, 3, ... n} is the set of customer points, K = {1, 2, ... k} is the set of available vehicles in the warehouse, i, j represents a single warehouse or customer point, i, j ∈ N, i ≠ j; represents the unit distance cost of vehicle k traveling from customer point i to customer point j, d ij represents the actual distance between node i and node j, α1 represents the penalty cost of early vehicle arrival, and α2 represents the penalty cost of late vehicle arrival; e i represents the left time window of customer i, t i represents the time when the vehicle arrives at customer i, l i represents the right time window of customer i; is a decision variable 0-1 variable, which takes 1 when vehicle k is traveling on path (i, j), otherwise it takes 0

[0100] S82: The objective function setting constraints of the pharmaceutical logistics vehicle routing optimization model with time window constraints described in S8 are as follows:

[0101]

[0102] τ0=e0 (6)

[0103]

[0104] Formula (1) indicates that the cold storage center has enough vehicles to complete the delivery task; Formula (2) indicates that every customer is served; Formula (3) indicates that each vehicle starts from the cold storage center and finally returns to the cold storage center; Formula (4) indicates that the total demand of customers served by each vehicle must not exceed the maximum load of the vehicle; Formula (5) indicates the demand of customer point i is met; Formula (6) indicates that the departure time of the vehicle is the earliest opening time of the distribution center; Formula (7) indicates that the time when the vehicle returns to the distribution center must be earlier than the closing time of the distribution center; Formula (8) indicates that the time when the vehicle reaches the next node is equal to the sum of the time when it arrives at the previous node, the service time of the previous node and the travel time between the two nodes; Formula (9) is the value constraint of the decision variable.

[0105] S9: Design a genetic algorithm improved by the adaptive large-area search algorithm. In the first stage, the greedy algorithm is used to generate the initial solution. In the second stage, the genetic algorithm is used to optimize the generated initial solution. In the third stage, the adaptive large-area search algorithm is used for deep optimization. The solution optimized by the genetic algorithm is used to avoid falling into the local optimum of the genetic algorithm to obtain the optimal distribution plan. GA-ALNS (genetic algorithm improved by the adaptive large-area search algorithm) is used to solve the objective function of the vehicle path optimization of pharmaceutical logistics with time window constraints. The specific steps are as follows:

[0106] S91: Use the greedy algorithm to generate the initial solution: select a customer point i that is closest to the distribution center in the customer set, add it to the individual, and then search among the customers who have not been added to the individual to find the customer j that is closest to the current customer and add it to the individual, and repeat this operation until all customers are added to the individual.

[0107] Since the transport vehicle has a weight limit, the algorithm should follow the following principles when generating the initial population: when the total demand of customers is less than the maximum load of the vehicle when the number of customers is accumulated to n+1, the result is the opposite. In this case, the number of codes of the distribution center is inserted after the nth customer, indicating that the coding of a chromosome is completed. For example, if there are 5 customers in total, after the first and second customers are inserted, the third customer is inserted, which violates the weight constraint. In this case, the chromosome code is 1-2-6, and this is repeated until all individuals are included to form a number that can meet the population size.

[0108] S92: Calculate the fitness value of the population. The higher the fitness value of an individual in the population, the better the performance of the individual. Since the objective function set in this paper is to find the minimum value, the fitness function is taken as the inverse of the objective function.

[0109] S93: Perform selection, crossover, and mutation operations on all individuals in the population to obtain a new population P1. The selection operation here is a binary tournament selection operation, the crossover operation is an OX crossover method, and the mutation operation is a two-point mutation.

[0110] S94: Randomly select R chromosomes from the new population P1 as another population P1.1, and the remaining chromosomes are recorded as P1.2.

[0111] S95: The population P1.1 is recorded as Li and enters the ALNS operation. First, the ALNS parameters are initialized, and then all chromosomes in Li are destroyed and reorganized through the destruction operator and repair operator in ALNS.

[0112] Remove the operator set h r as follows:

[0113] (1) Similarity removal operator: remove the node with the smallest correlation with the selected node, that is, randomly select a customer i to remove, calculate the correlation D between the remaining customers and i, and select the customer with the largest correlation to remove, that is,

[0114]

[0115] Where V ij Is i and j on the same road, that is, whether they are served by the same vehicle. If i and j are on the same road, it is 1, otherwise it is 0; d inRepresents the distance from customer i to any customer.

[0116] (2) Random loop removal operator: selects a path with the largest total cost and then randomly selects m customers in the path to remove.

[0117] (3) Worst distance removal operator: Calculate the sum of the distances between each customer and the customers before and after him on the vehicle path, and select the top m customers with the largest sum of distances.

[0118] (4) Worst time removal operator: select the customer with the largest time deviation and remove it from the current path. The time deviation of customer j WT(j) = |y j -a j |, where y j is the time when the delivery vehicle arrives at customer j, a j is the earliest acceptable service time for customer j.

[0119] Insert operator set h i as follows.

[0120] (1) Optimal greedy insertion operator: Find feasible insertion locations for the destroyed customer points one by one, that is, satisfy the time window and load constraints, and select the location with the smallest path cost increment among all feasible insertion locations for insertion.

[0121] (2) Random insertion operator: randomly insert any position of the solution after destruction under the premise of satisfying the time window and load constraints.

[0122] (3) Optimal distance insertion: Find feasible insertion locations for the destroyed customer points one by one, that is, satisfy the time window and load constraints, and select the location with the smallest distance increment from all feasible insertion locations for insertion.

[0123] (4) Best time insertion heuristic operator: For each subpath in the solution to be reconstructed, under the constraints of time window and capacity, first find all feasible insertion positions for customer j and calculate the time deviation after insertion. The time deviation of customer j, BT(j) = |y j -a j |, where: y j is the time when the delivery vehicle arrives at customer j, a j is the earliest acceptable service time for customer j, and then customer j is inserted at the position with the smallest time deviation.

[0124] S96: When performing repair and destruction, each time a roulette wheel is used to select one operator from the four destruction operators and four repair operators for destruction and repair, therefore, an adaptive mechanism is designed to evaluate the effects of operators by scoring and change the probability of selecting them. The specific operation is as follows:

[0125] The weight of the current operator participating in the roulette selection is calculated based on the average score of each round. Initially, the weight of each operator is set to 1, and the initial score is set to 5. If the solution found by the operator is better than the current global optimal solution, the operator will add α1 points; if the solution found by the operator is better than the current solution, the operator will add α2 points; if the solution found by the operator is worse than the current solution, the operator will add α3 points. In addition, a weight update coefficient ρ is set to avoid falling into the local optimum due to too fast convergence speed. The weight of the operator is updated according to the following formula:

[0126]

[0127] Where: a is the weight, s a is the cumulative score, μ a The cumulative number of times used.

[0128] S97: Use the following rules to control the update of the solution after ALNS. Let s1 represent the current solution, s2 represent the new solution, f(s1) and f(s2) be the objective function values ​​of s1 and s2 respectively. If f(s2)-f(s1)<0, replace the current solution with the new solution, otherwise reject the current solution. ALNS ends when all chromosomes in Li have been ALNSed.

[0129] S98: Reorganize P1.1 and Li after ALNS operation into a population P1. Because in the selection operation, one duplicate individual is retained and the rest are deleted. Therefore, after crossover, mutation, and ALNS operations, the number of new offspring is less than the original population. It is necessary to use recombination operation to combine P1 with the initial population to keep the population unchanged.

[0130] S99: Determine whether the maximum number of iterations has been reached. If not, return to the second step until the maximum number of iterations has been reached.

[0131] Finally, it should be noted that the above embodiments are only used to illustrate the technical solution of the present invention rather than to limit it. Although the present invention has been described in detail with reference to the preferred embodiments, those skilled in the art should understand that the technical solution of the present invention can be modified or replaced by equivalents without departing from the purpose and scope of the technical solution, which should be included in the scope of the claims of the present invention.

Claims

1. The pharmaceutical logistics vehicle path optimization method based on a two-stage heuristic algorithm according to claim 1 is characterized by: The method consists of two phases; Phase 1: Use the map API to collect the longitude and latitude information of a certain route of a pharmaceutical company, including warehouses and all customer points, and convert the longitude and latitude information into metric coordinates; Randomly select K points as the initial cluster centers, where in represents the total customer demand, Q max represents the maximum load of the vehicle, i represents the customer point, and n represents the number of customers; Calculate the spatiotemporal distance between the remaining customer locations and the cluster center, where the spatiotemporal distance includes spatial distance and temporal distance. The spatial distance is calculated by Euclidean distance, and the temporal distance is calculated by the satisfaction degree of the service time window. A penalty coefficient is introduced for the temporal distance that exceeds the time window range. According to the vehicle load constraint, the remaining customers are assigned to the partition corresponding to the nearest cluster center according to the time and space distance; If a customer does not meet the current load constraint, it will be assigned to another partition with the smallest time and space distance from it; Check whether there are any customers who have not been delivered. If so, jump to the step of "based on the vehicle load constraint, assign the remaining customers to the partition corresponding to the nearest cluster center according to the time and space distance". Otherwise, recalculate the cluster center and repeat the above steps until the clustering results converge. Phase 2: Use the map API to calculate the distance matrix between the warehouse and all customer points, and export the distance matrix as a CSV file; Construct a pharmaceutical logistics vehicle routing optimization model, with the optimization goal of minimizing transportation cost and time window penalty cost, and setting vehicle capacity limit and time window constraint conditions; Design a genetic algorithm combined with an improved adaptive large domain search algorithm, including: Generate an initial solution using a greedy algorithm; Calculate the fitness value of individuals in the population; Perform selection, crossover and mutation operations on the population to obtain a new population; Randomly select some chromosomes from the new population to perform adaptive large-area search algorithm operations, including: Initialize ALNS parameters; Destroy and reorganize the chromosomes, where the destruction operators include similarity removal operator, random loop removal operator, worst distance removal operator, worst time removal operator, and the repair operators include optimal greedy insertion operator, random insertion operator, best distance insertion operator, and best time insertion heuristic operator; Adaptively adjust according to the operator effect to control the update of the solution; The chromosomes that have undergone ALNS operation and the remaining chromosomes are reorganized into a new population; Determine whether the maximum number of iterations has been reached. If not, return to calculate the fitness values ​​of individuals in the population until the maximum number of iterations is reached.

2. The pharmaceutical logistics vehicle path optimization method based on a two-stage heuristic algorithm according to claim 1 is characterized by: In the adaptive large domain search algorithm, the weight of the operator is updated according to the following formula: Where: a is the weight, s a is the cumulative score, μ a is the cumulative number of times used, and ρ is the weight update coefficient in the adaptive mechanism.

3. The pharmaceutical logistics vehicle path optimization method based on a two-stage heuristic algorithm according to claim 1 is characterized in that: In the adaptive large domain search algorithm, if the objective function value of the new solution is smaller than the objective function value of the current solution, the new solution replaces the current solution.

4. The pharmaceutical logistics vehicle path optimization method based on a two-stage heuristic algorithm according to claim 1 is characterized by: In the adaptive large domain search algorithm, if the objective function value of the new solution is greater than or equal to the objective function value of the current solution, the new solution is rejected.

5. The pharmaceutical logistics vehicle path optimization method based on a two-stage heuristic algorithm according to claim 1 is characterized by: When the greedy algorithm generates the initial solution, when the total customer demand is accumulated to the nth customer, the total customer demand is less than the maximum vehicle load. When the total customer demand is accumulated to n+1, the result is the opposite. At this time, the code number of the distribution center is inserted after the nth customer, indicating that the coding of one chromosome is completed.

6. The pharmaceutical logistics vehicle path optimization method based on a two-stage heuristic algorithm according to claim 1 is characterized by: In the genetic algorithm, the selection operation is a binary tournament selection operation, the crossover operation is an OX crossover method, and the mutation operation is a two-point mutation.

Images