Designing social distancing policies for the COVID-19 pandemic: A probabilistic model predictive control approach

1.   Introduction

The logistics industry has strongly promoted the development of the national economy under the impetus of information technology. However, this has led to high energy consumption and carbon emissions. Only by insisting on economic development, clean development and safe development, that is, sustainable development, can we achieve sound and rapid economic development. Sustainable development is an economic model that focuses on long-term development, while energy conservation and emission reduction are necessary measures for sustainable development, and the two complement each other.

Energy consumption of the logistics industry in China accounts for 10% of the total national energy consumption ^[1], which is not conducive to the development concept of green GDP ^[2]. Green GDP is the cost of economic loss caused by environmental pollution, degradation of natural resources, poor education, uncontrolled population, and poor management of GDP ^[3]. Therefore, the logistics industry is facing the formidable pressure of energy saving and emission reduction. The logistics system established with the help of modern information technology can effectively improve operational efficiency, and achieve cost reduction, energy saving, and consumption reduction ^[4]. Specifically, multi-center joint distribution can break the limitation of the distribution area, reduce the round trip, improve the efficiency of distribution, and achieve the goal of energy saving and emission reduction, which is an important choice for the low-carbon and green development of logistics enterprises ^[5].

Malik and Kim ^[6] considered the production rate affects carbon emissions generation in production, i.e., generally, higher production rates result in more emissions. Sarkar et al. Introduced a three-echelon sustainable supply chain model with a single supplier, single manufacturer, and multiple retailers. In this supply chain, the main task of the manufacturer is to reduce defective products and control carbon emissions to maintain sustainability. Numerical experiments find that the model obtains the global optimum solution at the optimum values of the decision variables and the reduction of carbon emission has been proved ^[7].

Ghanbarpour and Gustafsson ^[8] proposed that customers' positive perceptions of firm actions do not directly impact financial earnings; however, they do impact earnings through customer satisfaction. Tayyab et al. analyze a cleaner multi-stage production management system for carbon emissions reduction and active participation in corporate social responsibility activities, while also advancing the system economically. According to the model results, the lower batch size is the optimal strategy in response to increasing holding and CSR activity costs, but otherwise, for the setup cost in a serial production system ^[9].

Dantzig and Ramser ^[10] proposed the vehicle routing problem (VRP) in 1959. Since its inception, VRP has sparked considerable interest. In particular, the research findings on VRP considering carbon emission, time window and customer satisfaction have been fruitful. In order to solve minimum total travel time and minimum fuel consumption in VRP, Cinar et al. used a variable travel distance and load weight to construct the cumulative function, and used K-means clustering and C & W algorithm to form a two-stage algorithm to improve the algorithm search capability ^[11]. Xiao et al. considered vehicle type, speed, load weight, and time window in the vehicle path optimization process, and the optimization can reduce carbon emissions by 8% ^[12]. Deng et al. considered both the soft time windows and hard time windows and established a multi-modal transport route optimization model with the objectives of minimizing total cost and carbon emissions ^[13]. Zhang et al. introduced a cold chain logistics route optimization model considering the low-carbon economy, an ant colony optimization was used to avoid the influence of unreasonable parameter selection of algorithm performance ^[14]. Tao et al. proposed a distribution route optimization method to minimize carbon emission cost as well as comprehensive cost ^[15]. Zhang et al. proposed a multi-objective route optimization model of the instant distribution system taking into account total cost and customer time satisfaction, testing results showed that a slight rise in delivery cost can achieve more on-time delivery demand ^[16]. Li and Zhou presented a logistics distribution center location model which considered the influence of customer satisfaction, dynamic and static carbon emission, construction and operation cost, and proposed penalty cost under a time window to measure customer satisfaction ^[17].

With the development of green logistics, a single distribution center can no longer meet the distribution needs, and scholars have gradually conducted research on the joint distribution model of multiple distribution centers. Rabbani et al. presented a deterministic vehicle routing problem model taking into account multiple middle depots ^[18]. Adelzadeh et al. studied the vehicle routing problem of multi-distribution center joint distribution under time window constraints, providing a reference value for enterprises to reduce their distribution costs and improve their service level ^[19]. Li et al. established an optimization model of joint distribution terminal distribution path with carbon emission, taking into account various costs with a time window, with the total distribution cost optimized. The model was solved using a genetic algorithm ^[20]. Golestani et al. investigated a green hub location problem with the objectives of minimizing the system's total cost and maximizing the quality of the delivered product to the customer ^[21].

At present, there exists a number of literature on the research of carbon emission optimization, analyzing the role of low carbon and other factors in the distribution systems, and using intelligent optimization algorithms such as heuristic algorithms and exact algorithms ^[22,23,24] to solve similar VRP problems ^[25,26], and some also study the customer satisfaction problems in logistics and distribution processes such as centralized distribution, joint distribution, and cold chain distribution. However, there are few logistics distribution models that consider both carbon emissions and customer satisfaction. Logistics enterprises, can reduce carbon emissions through vehicle path optimization and achieve cost reduction. At the same time, economic benefits play a critical role. Therefore, it is meaningful to achieve both economic and environmental benefits for enterprises by maximizing customer satisfaction and reducing carbon emissions in the distribution process. Based on the above analysis, this paper presents a joint distribution network with multiple distribution centers, constructing a multi-objective planning model taking into account distribution cost, customer dissatisfaction, and carbon emission. An improved ant colony optimization (ACO) is proposed to solve the model by using an elite strategy and saving strategy to obtain an optimal solution.

2.   Multi-center joint distribution optimization model

2.1.   Problem description

The logistics industry has a long history of high energy consumption and carbon emissions ^[27]. Reducing energy consumption and carbon emission in the logistics process helps enterprises to realize both economic and social benefits. This paper solves the problem of joint distribution with multiple distribution centers and customer time window requirements, i.e., it is an MDVRPTW problem. In the multi-centers joint distribution network, the customer order is no longer handled by a fixed distribution route, but the appropriate distribution vehicle and distribution route are decided according to the overall planning of the logistics system. This distribution method can provide rapid response to customer time window demand and maximize delivery efficiency. The sharing of vehicles and inventory resources can effectively solve the problems of urban traffic congestion and vehicles running empty. The vehicle selects the original distribution center for replenishment, then continues to deliver until all orders are completed, and finally returns to the original distribution center. In this paper, carbon emissions and customer dissatisfaction are firstly measured. The carbon emission in operation is calculated by the vehicle's actual load rate and vehicle driving distance, and the customer dissatisfaction is calculated by the fuzzy mathematical affiliation function. Secondly, the joint distribution model of multiple distribution centers considering carbon emission and time window is constructed. The model takes into account actual vehicle loading rate, vehicle capacity, driving speed, loading and unloading time, number of vehicles, and number of distribution centers in the distribution network. As a result, this paper presents a multi-objective path optimization model with the objectives of minimum cost, minimum customer dissatisfaction, and minimum carbon emission.

Assumptions in this work are listed as follows.

1) The location of the distribution center is known, the reserve capacity can well serve each demand point, and the transport capacity is able to complete the distribution task.

2) Each vehicle departs from its respective distribution center and returns to the same departure point after completing the distribution task.

3) The vehicle type is consistent.

4) Each node is served by one vehicle once only, and the node requirements can be satisfied in a one-time service.

5) The time window requirements of the node must be satisfied.

6) The loading quantity of a vehicle can satisfy the total demand in its corresponding path.

7) Each vehicle is not reassigned in the middle of the distribution.

2.2.   Parameter description

The meaning of the parameters in the model are presented as follows.

M: vehicle set in distribution centers.

${V}_{m}$: orders served by vehicle m.

$N$: order set.$\left|N\right| = n$, there are $n$ orders.

$i, j$: indicates order number, also used to mark customer number of the order.

${d}_{ij}$: distance from customer$i$ to $j$, ${d}_{ij} = {d}_{ji}$.

${x}_{ij}^{m}$: decision variable, whether vehicle m serves customer $j$ directly after serving customer $i$.

${y}_{i}^{m}$: decision variable, whether customer $i$ is served by vehicle m or not.

$Q$: vehicle loading capacity and capacity of all vehicles are the same.

${d}_{i}$: demand of customer $i.{d}_{0} = 0$, which means the demand of the distribution center is 0.

TD: the maximum traveling distance of the vehicle.

${t}_{i}$: the arrival time of vehicle m reaches customer $i.$

${\mu }_{i}\left({t}_{i}\right)$: dissatisfaction rate of customer $i$ when served at time ${t}_{i}$.

${p}^{fuel}$: unit fuel price.

${c}_{p}^{e}$: penalty cost per unit of time if the vehicle arrives early.

${c}_{p}^{l}$: penalty cost per unit of time if the vehicle arrives late.

${C}_{1}$: distribution cost.

${C}_{2}$: customer dissatisfaction.

${C}_{3}$: carbon emission.

2.3.   Distribution cost calculation

The costs generated by distribution are composed of vehicle driving fuel cost, and penalty cost if the vehicle arrives early or late.

${\mathrm{C}}_{1}^{1}$ presents vehicle driving fuel cost, which is calculated by following equation.

The penalty cost consists of early and late arrival penalty cost. The early arrival penalty cost is ${\mathrm{T}}_{\mathrm{m}, \mathrm{i}}^{\mathrm{e}} = \mathrm{m}\mathrm{a}\mathrm{x}\left(0, {\mathrm{a}}_{\mathrm{i}}-{\mathrm{t}}_{\mathrm{i}}\right)\hspace{0.33em}, \hspace{0.33em}\mathrm{m}\in \mathrm{M}, \hspace{0.33em}\mathrm{i}\in \mathrm{V}$ and the late one is ${\mathrm{T}}_{\mathrm{m}, \mathrm{i}}^{\mathrm{l}} = \mathrm{m}\mathrm{a}\mathrm{x}\left(0, {\mathrm{t}}_{\mathrm{i}}-{\mathrm{b}}_{\mathrm{i}}\right)\hspace{0.33em}, \hspace{0.33em}\mathrm{m}\in \mathrm{M}, \hspace{0.33em}\mathrm{i}\in \mathrm{V}$. Therefore, the penalty cost can be calculated as following equation.

At last, the vehicle distribution cost can be obtained from formula (3).

2.4.   Customer dissatisfaction calculation

The VRPTW only considers how to complete distribution tasks with the least number of vehicles and the shortest driving distance to save cost. However, in the retail industry, with the diversification of cargo and the increase of real-time distribution requirements from customers, multiple small-batch distributions are becoming more and more popular. Customers make delivery appointments at any time, and they have a certain tolerance for early or late arrival. In the actual situation, untimely delivery, in the long run, will damage enterprises' reputations. Consequently, considering customer dissatisfaction with delivery efficiency is of great practical importance.

The method of determining customer dissatisfaction is shown in Figure 1. ${\mathrm{t}}_{\mathrm{i}}$ is the vehicle arrival time at customer $\mathrm{i}$. The agreed time window is $\left[{\mathrm{ET}}_{\mathrm{i}}, \mathrm{L}{\mathrm{T}}_{\mathrm{i}}\right]$. If the vehicle arrives at the customer's during this period, customer dissatisfaction is 0. The tolerable early arrival time of the customer is $\mathrm{E}{\mathrm{T}}_{\mathrm{i}}-{\mathrm{a}}_{\mathrm{i}}$, and the tolerable delayed arrival time is${\mathrm{b}}_{\mathrm{i}}-\mathrm{L}{\mathrm{T}}_{\mathrm{i}}$. If the vehicle arrives at the customers within the tolerable time window, although the order still would be received, the customer would be dissatisfied depending on the degree of order advancement or delay. The customer's intolerable early arrival time is${\mathrm{a}}_{\mathrm{i}}$, and the intolerable delay time is ${\mathrm{b}}_{\mathrm{i}}$. If the vehicle arrives outside of the time period $[{\mathrm{a}}_{\mathrm{i}}, {\mathrm{b}}_{\mathrm{i}}]$, customer dissatisfaction is 1.

Therefore, based on above analysis, when customer$\mathrm{i}$receives order at time ${\mathrm{t}}_{\mathrm{i}}$, the dissatisfaction rate ${\mathtt{μ }}_{\mathrm{i}}\left({\mathrm{t}}_{\mathrm{i}}\right)$ can be presented as follows.

Figure 2 shows customer dissatisfaction when receiving orders at different moments.

Customer dissatisfaction consists of the dissatisfaction of all customers, i.e., total customer dissatisfaction rate is

2.5.   Carbon emission calculation

2.5.1.   Fuel consumption rate (FCR) calculation

FCR during vehicle driving is closely related to the actual load factor of the vehicle, this factor is an important indicator to measure the effective utilization of the vehicle. The higher the actual load factor of the vehicle, the higher the cost of fuel consumption, but the lower the transportation cost, the higher the utilization rate of the vehicle. According to the FCR regression equation ^[28], the fuel consumption per unit distance can be calculated as below.

Where, $\mathtt{ε }$ is the fuel consumption per unit distance, $\mathrm{Q}$ is vehicle load capacity, and ${\mathrm{Q}}^{\mathrm{*}}$ is the maximum vehicle load. When the vehicle is empty, ${\mathtt{ε }}_{0} = 0.254$, when it's fully loaded, ${\mathtt{ε }}^{\mathrm{*}} = 0.37944$.

Let actual load factor ${\mathrm{p}}_{\mathrm{i}} = \frac{\mathrm{Q}}{{\mathrm{Q}}^{\mathrm{*}}}$, FCR is displayed in Eq (7) ^[29].

2.5.2.   Vehicle fuel volume calculation

${\mathrm{d}}_{\mathrm{i}\mathrm{j}}$ is the traveling distance between customer $\mathrm{i}$ to customer $\mathrm{j}$, ${\mathrm{Q}}_{\mathrm{i}\mathrm{j}}$ is the vehicle load during this driving interval, $\mathtt{ε }\left({\mathrm{Q}}_{\mathrm{i}\mathrm{j}}\right)$ is FCR, and ${\mathrm{c}}_{\mathrm{i}\mathrm{j}}^{\mathrm{f}\mathrm{u}\mathrm{e}\mathrm{l}}$ is the fuel consumption, then ${\mathrm{c}}_{\mathrm{i}\mathrm{j}}^{\mathrm{f}\mathrm{u}\mathrm{e}\mathrm{l}}$ can be expressed as below ^[30].

2.5.3.   Carbon emission calculation

Emissions of CO₂ during vehicle operation and fuel consumption ${\mathrm{c}}_{\mathrm{i}\mathrm{j}}^{\mathrm{f}\mathrm{u}\mathrm{e}\mathrm{l}}$ on the combustion conversion factor $\mathrm{F}$ is linearly related ^[31]. Carbon emission is measured by Eq (9).

Different fuels have different combustion conversion factors. This work assumes that minivans are used for distribution, according to literature ^[32], the combustion conversion factor F = 2.3, which means that the weight of CO₂ released per liter of gasoline combustion is 2.3 kg.

2.6.   Model formulation

Based on the above sub-item analysis, the multi-objective optimization model for joint distribution considering carbon emission and customer dissatisfaction is constructed by using a dictionary order multi-objective planning method.

Equation (10) is the multi-objectives.

Equation (11) is the priority of different objectives, which can be adjusted according to the objective priority during decision-making.

The loading capacity of the vehicle and the maximum driving distance are considered during distribution.

s.t.

Constraint (12) indicates that vehicle m can only serve one customer $\mathrm{i}$ before serving customer $\mathrm{j}$. Constraint (13) means that after serving customer $\mathrm{i}$, vehicle m should serve customer $\mathrm{j}$ immediately. Constraint (14) indicates that each customer can only be served by one vehicle. Constraint (15) indicates that the total distance traveled by the vehicle does not exceed its maximum traveling limit. Constraint (16) presents that the sum of customer demands delivered by the vehicle does not exceed its maximum capacity. Constraint (17) is to avoid loops during delivery. Constraint (18) means the time constraint of the vehicle at each customer point.${\mathrm{x}}_{\mathrm{i}\mathrm{j}}^{\mathrm{m}}$ and ${\mathrm{y}}_{\mathrm{i}}^{\mathrm{m}}$ are decision variables. When vehicle m serves customer $\mathrm{j}$ immediately after serving customer $\mathrm{i}$, ${\mathrm{x}}_{\mathrm{i}\mathrm{j}}^{\mathrm{m}} = 1$, otherwise, ${\mathrm{x}}_{\mathrm{i}\mathrm{j}}^{\mathrm{m}} = 0$. When customer $\mathrm{i}$ is served by vehicle m, ${\mathrm{y}}_{\mathrm{i}}^{\mathrm{m}} = 1$, otherwise, ${\mathrm{y}}_{\mathrm{i}}^{\mathrm{m}} = 0$.

3.   Improved ant colony algorithm

This paper studies a multi-objective optimization problem, which is an NP-HARD problem, and heuristic algorithms are widely used in solving such problems ^[33]. Among them, the ant colony algorithm uses a distributed parallel computing mechanism, which has the advantages of positive feedback, good robustness, and obtaining a more satisfactory and feasible solution within an acceptable time frame. The basic idea of the ant colony algorithm: ants release pheromones along the way during crawling, and the pheromone will gradually weaken with time. When other ants pass the path, the pheromone is strengthened, and the ants tend to crawl in the direction of high pheromone intensity during the crawling process ^[34]. The more ants walk on a certain route, the higher the concentration of pheromone left on the route. A large number of ants form a positive information feedback mechanism, and finally, a route with the largest number of crawling ants is formed, which is also the optimal route ^[35].

However, it has disadvantages of early maturity and low solution efficiency. Therefore, an improved ant colony algorithm is proposed in this paper.

3.1.   Coding scheme

Considering both vehicle load and farthest driving distance in joint distribution is a new VRP. The taboo table design is different from the traditional solution. It needs to find the original warehouse for replenishment, and then continue. What's more, the vehicle has a limit of maximum traveling distance. In order to solve the problem described above, this paper designed a specific taboo table coding method of ACO, which is shown in Figure 3.

For example, in Figure 3, numbers 23, 24 and 25 represent warehouses, and numbers 1 to 22 represent customers.

The route of a certain ant contains two distribution vehicles, where the driving path of vehicle 1 is 23→2→11→15→5→7→23→13→17→1→23→18→21→9→10→16→23, while the driving path of vehicle 2 is 23→6→12→20→19→23→3→8→14→4→22→23. Distribution vehicle 1 departs from warehouse number 23 and delivers orders according to route 23→2→11→15→5→7. After completing the delivery task to customer 7, there are no more goods for the next customer, therefore, vehicle 1 goes to the original warehouse number 23 for replenishment, after reloading, it continues delivery along route 23→13→17→1. And goods are out again at customer 1, so vehicle 1 returns to the original warehouse number 23 for reloading, and then continues distribution until the completion of all orders assigned to vehicle 1. Vehicle 2 adopts a similar distribution method as vehicle 1. The driving routes of all vehicles constitute an ant's route, and the routes of all ants constitute a taboo table.

3.2.   Vehicle selection rules

The distribution center needs to select a vehicle for delivery orders. Since the demand of each customer is unequal, and the distance between two adjacent customers is different, in order to avoid uneven free and busy times of vehicles, the distribution center gives priority to the vehicle with lower busyness for delivery according to a certain probability.

Because it is in accordance with the probability to choose the distribution vehicle, there will be some vehicles arranged for more distribution tasks, and some less. For example, vehicle 1 delivers more goods to customers far away, while vehicle 2 delivers fewer goods to customers nearer, then, vehicle 1 takes more time and vehicle 2 takes less. In order to avoid this problem, the busy level of distribution vehicles is determined in accordance with the current distribution time of each vehicle. Such as, if vehicle 1 spends more time, and vehicle 2 spends less, then vehicle 1 is busier than vehicle 2, and vehicle 2 will be given priority for new assignments. The vehicle selection rules are shown in Figure 4.

3.3.   Customer selection rules

After completing vehicle selection, it is time to assign customers to be served by specific vehicles under constraints of maximum load capacity, the farthest driving distance, and other factors. When assigning customers, the transfer probability of the vehicle traveling from the current customer to the remaining customers is firstly calculated according to the indexes such as pheromone and visibility, then the customer with maximum transfer probability is selected. If the customer with the maximum transfer probability is not selected, then roulette is constructed according to the transfer probability to select the customer to be served. After completing customer selection, whether the customer's demand exceeds the maximum load capacity of the vehicle is judged, if it exceeds the maximum load capacity, then go back to the last step, if it still does not meet the maximum load capacity of the vehicle, then select the original warehouse to the current customer for replenishment, and continue to select the customer to visit according to the probability after replenishment. Customer selection rules are shown in Figure 5.

3.4.   Saving strategy

In the transfer rule of the basic ant colony algorithm, the transfer probability takes into account the pheromone concentration and visibility factors. After departure, whether a distribution vehicle visits point $\mathrm{j}$ from point $\mathrm{i}$ is determined by the probability. In order to accelerate the convergence speed of the ant colony algorithm, the maximum transfer probability customer point is selected with a certain probability, and if the probability is not satisfied, the next point to be visited is selected according to the probability of each point. The pseudo-random rule is: rand is a random number between 0 and 1, when$\mathrm{rand} < = \mathrm{m}\mathrm{a}\mathrm{x}\mathrm{s}\mathrm{h}\mathrm{i}\mathrm{f}\mathrm{t}\_\mathrm{p}$, the next point to be visited is selected according to the maximum probability, otherwise a probability roulette is constructed and the point to be visited is selected according to the roulette. The probability of transfer from point $\mathrm{i}$ to point $\mathrm{j}$ is:

In order to improve the searchability of the basic ant colony algorithm, the saving factor is added to the transfer probability by using the idea of the saving algorithm as a reference, and the transfer probability is presented in Eq (22).

$\mathrm{a}\mathrm{l}\mathrm{l}\mathrm{o}\mathrm{w}\mathrm{e}{\mathrm{d}}_{\mathrm{m}}$ indicates the set of points that m vehicle can be accessed from current point $\mathrm{i}$. Still, $\mathtt{α }, \mathrm{ }\mathtt{β }$ indicate the importance of pheromone, and the importance of visibility respectively, and ${\mathtt{τ }}_{\mathrm{i}\mathrm{j}}$ indicates the pheromone concentration from point $\mathrm{i}$ to point $\mathrm{j}$, ${\mathtt{η }}_{\mathrm{i}\mathrm{j}}$ is the visibility from point $\mathrm{i}$ to point $\mathrm{j}$. ${\mathtt{μ }}_{\mathrm{i}\mathrm{j}}$ is the amount of savings of edge$\left(\mathrm{i}, \mathrm{j}\right)$, i.e., ${\mathtt{μ }}_{\mathrm{i}\mathrm{j}} = {\mathrm{h}}_{\mathrm{i}\mathrm{b}}+{\mathrm{h}}_{\mathrm{b}\mathrm{j}}-{\mathrm{h}}_{\mathrm{i}\mathrm{j}}\left(\mathrm{b}\in \mathrm{D}\right)$. The savings of the edge reflects the importance of the path$\left(\mathrm{i}, \mathrm{j}\right)$ ^[36].

3.5.   Pheromone update

After the ant finishes constructing paths of all delivery vehicles, it needs to update the pheromones on the paths. Since the problem is a multi-objective optimization with cost as the first priority, and the cost is mainly affected by driving distance and full load rate, therefore, the pheromone update for the ant $\mathrm{k}$ should refer to travel distance ${\mathrm{L}}^{\mathrm{k}}$ of ant $\mathrm{k}$. The pheromone update equation is as the following equation.

Where $\mathtt{ρ } \; \left(0 < \mathtt{ρ } < 1\right)$ is the pheromone volatility factor, $△{\mathtt{τ }}_{\mathrm{i}\mathrm{j}}^{\mathrm{k}}$ is the pheromone increment of ant $\mathrm{k}$ from point $\mathrm{i}$ to point $\mathrm{j}$, $\mathrm{Q}$ is a constant, and ${\mathrm{f}}^{\mathrm{k}}$ is the multi-objective weighted sum of ant $\mathrm{k}$, and ${\mathtt{τ }}_{\mathrm{i}\mathrm{j}}$ is the pheromone concentration from point $\mathrm{i}$ to point $\mathrm{j}$.

The pheromone of the ant gradually evaporates during walking, and leaves a new pheromones, this process guides the optimal solution, which will be gradually strengthened and kept stable through repeated operations.

3.6.   Elite strategy

The ant colony algorithm relies on probability to guide the search when selecting the access point. It has been proved that the ant colony algorithm converges to the optimal solution with probability 1 after repeated iterations, and in order to speed up the convergence of the ant colony algorithm, the convergence of the algorithm can be accelerated by enhancing the pheromone of good ants ^[37]. In this paper, the elite individuals of each generation are preserved by the percentage of $\mathrm{e}\mathrm{l}\mathrm{i}\mathrm{t}\mathrm{i}\mathrm{s}\mathrm{t}\_\mathrm{c}\mathrm{o}\mathrm{e}\mathrm{f}$, which are merged with the next generation to enhance the number of excellent individuals and enhance the pheromone on the path of the excellent solution.

3.7.   Flow chart of improved ACO

The algorithm optimization process considers the busy degree of vehicles, gives priority to the less busy vehicles for distribution, and improves the local search ability of the algorithm through elite strategy. The vehicle selects the original warehouse for replenishment and uses the dictionary order method in multi-objective processing. The flow chart of the improved algorithm is shown in Figure 6.

4.   Numerical experiments

4.1.   Basic data

Using the simulation data, the validity of the model and the improved ACO are verified. In the set of this simulation experiment, enterprise M is a professional logistics company that provides service to stores and supermarkets. In this paper, 41 stores served by enterprise M are selected for data collection. The location coordinates of 41 stores and 3 warehouses, the time window requirements of customers, and the demand of each customer point are shown in Table 1.

As the operating hours of convenience stores are 8:00–22:00, in order to facilitate the management of goods, the ideal time to receive goods is 9:00–17:00. During the rest of the time, goods can be received at 8:00–9:00 and 17:00–22:00 with dissatisfaction generated due to staff shortage. Dissatisfaction of receiving goods before 8:00 and after 22:00 is 1. The time when the vehicle starts to deliver is set as 0, and the delivery time and time window are accumulated in turn.

4.2.   Parameter setting

The distribution process is set up in the distribution center with three distribution vehicles, each with a maximum load of 5 tons, and the proposed algorithm is implemented in MATLAB2017a. The parameters used in this paper are set according to the data of China Emission Accounts and Datasets and reference ^[38], as shown in the following Table 2.

4.3.   Performance analysis

Compare the improved ant colony algorithm and genetic algorithm with the basic ant colony algorithm. Both of them use joint distribution method, and there are 3 distribution vehicles. Take the best results of the three algorithms as the optimal solution for comparison. Table 3 displays the comparison of the number of replenishment, full load rate, cost, dissatisfaction and carbon emissions. As shown in Table 3, the improved ant colony algorithm has 95.49, 92.4 and 86.42% full load rates, which are better than the other two algorithms, respectively. It can also be seen that the cost and carbon emission of the improved ant colony algorithm is the smallest due to the highest vehicle loading efficiency of the improved ant colony algorithm, which also improves customer satisfaction.

Besides, the comparison reveals that the improved ant colony algorithm outperforms the basic ant colony algorithm and genetic algorithm in all metrics.

After the algorithm enters the stagnation phase, the first ranked ant route is selected as the satisfactory solution, and since the model deals with the actual problem, the satisfactory solution can be used as an effective solution for this problem. The solutions of this multi-objective optimization are shown in Table 4, while the carbon emission and customer satisfaction are shown in Table 5. According to the results shown in Table 4, the full load rate of all vehicles exceeds 85%, and the full load rate of vehicle 1 reaches 95.49%. The lowest cost obtained is 583.91 Yuan, which solves the cost loss problem to a greater extent; the lowest customer satisfaction is 1.948, and the vehicle full load rate improves the efficiency of cargo arrival, thus improving customer satisfaction; the lowest carbon emission is 45.15 Kg, and the increase of vehicle full load rate also leads to the reduction of the number of times the vehicle is used, thus suppressing the increase of carbon emission.

The routes of three delivery vehicles are obtained according to calculation results, which are shown in Figures 7–9, and the overall routes of three vehicles is presented in Figure 10.

4.4.   Discussion of results

From delivery routes presented in above figures, it is obvious to find that vehicle 1 is mainly replenished in distribution center 42 during the distribution process, with 2 replenishments. Vehicle 2 is mainly replenished in distribution center 43 with 1 replenishments. And vehicle 3 is mainly replenished in distribution center 44, with 3 replenishment. The cross-regional distribution characteristics of vehicles 1, 2 and 3 are clearly displayed in the distribution process.

The full load rate of the first vehicle is 95.49%, its carbon emission is 22.7031 KG, and the fuel cost is 59.23 Yuan. It has been replenished twice. The full load rate of the second vehicle is 92.4%, its carbon emission is 11.7817 KG, and its fuel cost is 30.73 Yuan. It is replenished once. The full load rate of the third vehicle is 86.42%, its carbon emission is 10.6698 KG, and its fuel cost is 27.83 Yuan. It has been replenished three times. Five of the 41 stores failed to deliver at the specified time, and their total dissatisfaction was 1.948.

5.   Managerial insights

This section discusses the industrial benefits for industry managers and are as follows:

The implementation of a reliable logistics and transportation system is a necessary condition for environmental protection. At the same time, a reliable logistics transportation system can also respond more quickly to the transportation route, reduce the transportation cost, and thus improve the core competitiveness of enterprises.

Based on above results analysis, it can be concluded that the solution of proposed algorithm achieves lower distribution cost with higher customer satisfaction when considering fuel consumption cost, penalty cost, and other factors in route planning. Considering the carbon emissions of trucks, it should pay attention to the full load ratio of trucks. The higher the carrying capacity, the more carbon emissions, the lower the carrying capacity, and the higher the cost. Consequently, VRP optimization with considerations of carbon emissions calculated by the actual weight of the vehicle and the driving mileage helps to realize green logistics.

For enterprises, multi-center joint distribution considering cost, customer dissatisfaction and carbon emissions can satisfy diversified customer needs, meeting requirement of multi-species, small-lot, and individualized orders and improve customer experience. This implementation makes contributions of reducing enterprise operation cost, improving logistics and distribution efficiency, as well as realizing good social value.

6.   Conclusions

To reduce damage to the environment during operation, logistics enterprises pay more attention to green logistics, thus, energy consumption and carbon emission as well as customer dissatisfaction are core issue to solve in distribution system. In this paper, a multi-objective planning model considering distribution cost, customer dissatisfaction and carbon emission is constructed, and an intelligent optimization algorithm is designed to solve the problem. Based on the characteristics of slow convergence and easy premature of traditional ant colony algorithm, elite strategy and saving strategy are used to improve the ant colony algorithm. The case study results show that the improved ACO has better performance compared to the basic ACO, which can effectively solve this multi-objective optimization problem and can be applied to similar multi-objective optimization. At the same time, the joint distribution of multiple distribution centers can break the barrier that resources cannot be shared between regions, realize the purpose of reducing costs and increasing efficiency of enterprises. In the process of multi-center joint distribution, considering distribution cost, customer dissatisfaction and carbon emissions can realize both economic and social benefits, which makes it be more in line with the requirements of green GDP regulations.

The future development direction of this paper is as follows:

Multi objective issues should be considered more comprehensively. This paper considers distribution cost, customer dissatisfaction and carbon emissions, but only considers distance, weight, full load ratio and other factors when quantifying these factors, lacking comprehensive consideration of more factors.

The solving efficiency of the improved ant colony algorithm should also be improved. The problem studied in this paper is NP hard problem, which is difficult to solve. When solving this multi-objective optimization problem, the ant colony algorithm is improved to improve the solution efficiency of the algorithm, but the search ability of the ant colony still needs to be enhanced.

Therefore, the next development direction is as follows:

The driving speed of the vehicle shall be considered when calculating carbon emissions.

When calculating customer dissatisfaction, consider the problem of traffic congestion when driving.

Vehicles are required to go to the nearest warehouse for replenishment.

The pheromone updating function in ant colony algorithm needs to be improved.

The proposed algorithm has no verification and exploration of multiple scenarios or multiple data sets.

Acknowledgments

This paper was supported by the Social Science Grand Research of the Hebei Education Department (ZD202105). The authors also gratefully acknowledge the helpful comments and suggestions of the reviewers, which have improved the presentation.

Conflict of interest

The authors declare that there are no conflicts of interest.

Data availability

The research data used to support the findings of this study are included within the article, the data used are shown in Tables 1 and 2.