# Resource Scheduling Method for Optimizing the Distribution Path of Fresh Agricultural Products under Low-Carbon Environmental Constraints

Authors: Qimin Fu; Jun Li; Huanhuan Chen
Journal: Scientific Programming
 (2022)

---

## Abstract

Aiming at the development concept of relative dispersion of urban multiregional business centers, saving energy, and emission reduction, this article proposes an algorithm for solving the problem of green logistics distribution to fresh agricultural products considering low-carbon and environmental protection, to realize the economic and environmental protection of cold chain logistics. Firstly, we comprehensively consider the fixed costs of delivery vehicles, transportation costs, damage costs of fresh agricultural products, cooling costs, carbon emissions costs, and penalty costs due to service time windows that do not meet customer requirements as the objective function to construct green logistics distribution location-routing optimization model. Then, Tent chaotic perturbation method is introduced to optimize the genetic algorithm and an improved genetic algorithm is used to solve the distribution location-routing optimization model to obtain the best scheme. Finally, the proposed algorithm is verified experimentally based on the MATLAB simulation platform. Experimental results show that the total distance and total cost of the distribution plan obtained by the proposed algorithm are 128.96 km and 12,593 yuan, respectively, which are superior to other comparison algorithms and can be used as a reference for logistics enterprises’ distribution decision-making.


---

## Body

## 1. Introduction

In recent years, greenhouse gas emission reduction in transportation has attracted much attention in the context of global climate change and environmental protection. In addition to vigorously developing new energy sources, the optimization of logistics and transportation routes is also an effective way to improve carbon emissions [1]. Especially for the cold chain logistics routing optimization with high energy consumption and high time efficiency requirements, it is of great significance for saving energy, reducing carbon emissions, and seeking sustainable economic and environmental development [2, 3]. Generally, a complete fresh cold chain logistics should consist of four links: refrigerated processing, refrigerated storage, refrigerated transportation and distribution, and refrigerated sales. Agricultural products should always be in a low temperature environment from the place of production to the final retail link, and the entire operation process is accompanied by the generation of energy consumption. Among them, refrigerated storage and refrigerated transportation are the links with the most energy consumption and the highest cost [4]. Therefore, optimizing energy saving and emission reduction in the transportation process will greatly promote the development of a low-carbon economy.The Vehicle Routing Problem (VRP) was first proposed by Dantzig and Ramser in 1959 and has since received extensive attention and in-depth research from all walks of life [5]. For general product logistics distribution routing optimization research, only the fixed cost and transportation cost of vehicle use need to be considered from the perspective of cost optimization. Usually, time window constraints are not considered, but the distribution cost of fresh and perishable food is different from the cost of logistics distribution at room temperature. Thus, the characteristics of cold chain logistics should be fully considered for research when optimizing the distribution path of cold chain logistics [6].Although the issue of the low-carbon cold chain has begun to receive attention in practice, there are not many academic research documents comprehensively focusing on cold chain and low-carbon topics. Some scholars qualitatively analyzed the challenges faced by cold chain logistics in a low-carbon environment and put forward some suggestions on energy saving and emission reduction from the perspective of the country and industry. However, in the research process of VRP, few scholars currently consider cold chain and low-carbon factors at the same time [7]. Therefore, in view of the green logistics planning problem, including low-carbon environmental protection, resource conservation, and time efficiency optimization, low-carbon environmental protection perspective is proposed to solve the problem of green logistics distribution location-routing optimization problem for fresh agricultural products under the perspective of low-carbon environmental protection. Compared with the traditional vehicle routing optimization algorithm, its innovations are as follows:(1)
Since most of the existing distribution routing optimization models consider economics and rarely involve carbon emissions, a green and low-carbon logistics distribution location-routing optimization model is constructed with the lowest cost as the objective function. It includes fixed costs, transportation costs, cargo damage costs, refrigeration costs, penalty costs, and carbon emission costs to achieve low-carbon distribution of fresh agricultural products and reduce energy consumption.(2)
In order to avoid the emergence of premature genetic algorithm and local optimum, the proposed algorithm introduces Tent chaotic perturbation method. This method uses a perturbation mechanism to initialize the population, performs perturbation again after the selection operation to increase the diversity of population, and designs the rules of selection, crossover, and mutation operators to speed up the solution.


## 2. Related Work

Due to the large scale of the city, the relatively scattered commercial centers, and the complex basic transportation, general urban logistics mostly adopts regional distribution. According to the size of the customer base, multiple distribution centers are set up to be responsible for the corresponding distribution services and so can abstract the optimal distribution problem into multidistribution center VRP (MDVRP) [8]. MDVRP is derived from VRP. After recent years of research, certain results have been achieved in the model establishment and solution methods, mainly focusing on precise algorithms and intelligent heuristic algorithms [9]. Among them, precise algorithms include branch and bound algorithms and dynamic programming algorithms. For example, [10] introduced an algorithm for generating linear independent feasible solutions for the problem of truck assignment and path planning at the docks to form an integer programming model. It determined the optimal allocation plan on the root node. Reference [11] proposed a framework model combining branching and cutting to optimize the planning of typical routes. This method effectively improved the efficiency of route planning and reduced the cost, but the operating cost is not considered for vehicle route planning. For example, [12] aimed at the dynamic path planning of unmanned environment monitoring vehicles under complex road conditions. Based on the idea of two-level programming, they proposed a hybrid algorithm that combines global and local path planning. It can effectively solve the path local optimization problem, but the optimization of path cost and carbon emission problems still needs to be in-depth. On the whole, there is relatively little research on accurate algorithm solving MDVRP model, and the solution effect is not satisfactory.Different from accurate algorithms, modern heuristic algorithms mainly use some intelligent algorithms to complete the routing optimization containing multiple goals, which is an important research method of MDVRP [13]. Currently, it mainly includes genetic algorithm, tabu search, quantum evolution, memetic algorithm, and deep learning algorithm. For example, [14] proposed a genetic algorithm-based motion trajectory planning method for industrial robots, aiming at the shortest time and lowest energy consumption to achieve the global optimal planning of the motion path. However, the planned path was biased towards ideal road conditions, the actual road conditions in the environment were more complicated, and the optimization capability needs to be improved. Reference [15] designed a search path planning algorithm based on an improved genetic algorithm to improve the adaptability of search path planning for the problem of Unmanned Aerial Vehicle (UAV) searching for targets with uncertain moving speed and direction. For example, [16] combined the advantages of unmanned aerial vehicles and proposed a new Voronoi path planning algorithm to solve the path planning problem of unmanned surface vehicles. Compared with the performance of the virtual machine algorithm, the proposed algorithm has better path planning performance. Reference [17] proposed a path planning model based on a convolutional neural network for the problem of unmanned vehicle path planning, which effectively improved the efficiency of vehicle path planning. However, there is still a lack of multiobjective optimization considerations for low-carbon and environmentally friendly green logistics path planning. Thus, a green logistics distribution location-routing optimization program from the perspective of low-carbon and environmental protection is proposed, which effectively solves the problem of location and routing optimization in the green and environmentally friendly distribution of fresh agricultural products.


## 3. Problem Description and Model Establishment

### 3.1. Problem Description

With the expansion scale of cities and the relative dispersion of commercial centers, there is an increasing demand for urban multiregional distribution and chain supermarket distribution. How to realize the effective integration of logistics resources and reduce the cost of urban logistics distribution has become an urgent problem to be solved in current urban logistics distribution [18]. The urban multiregional distribution problem is a typical MDVRP problem involving multiple distribution centers participating in the fulfillment of distribution requirements. The main difference between it and the routing optimization problem of a single distribution center is how different distribution centers in MDVRP jointly participate in the distribution process [19]. At present, large-scale vehicles in cities are affected by traffic restriction policy, and most of the enterprises’ logistics distribution adopts hierarchical distribution. Set up large warehouses on the outskirts of the city and deliver goods to various distribution centers through unrestricted small and medium-sized vehicles. The distribution centers adopt a partitioned distribution model to be responsible for customer service in each area, as shown in Figure 1.Figure 1
Schematic diagram of VRP in multiple distribution centers.MDVRP is generally described as follows: in the logistics network, multiple distribution centers are responsible for the distribution services of customers, and the locations of distribution centers and customers in the static MDVRP are known. At the same time, it is assumed that the goods in the warehouse and the distribution center can meet the needs of all customers, and the goods of different customer needs can be mixed with each other. It is required that the demand of each customer does not exceed the maximum load capacity of vehicles and that only one vehicle is responsible, and the demand is not allowed to be split [20, 21]. It is assumed that the maximum load of vehicles is known and overloading is not allowed. The maximum travel distance for one delivery is known, and violations are not allowed. The distribution distance between the distribution center and customers and between the customers is known. It is required that the vehicles start from the distribution center when the distribution task starts, and it is required to return to the original distribution center after the task is completed.


### 3.2. Model Establishment

#### 3.2.1. Fixed Cost of Vehicles

The fixed cost of using a vehicle is usually constant, including the fixed loss of vehicles and the cost related to the use of vehicles, such as the driver’s salary. It has nothing to do with the mileage of vehicles and the number of customers served. The fixed costΩ1 of vehicles is calculated as follows:(1)Ω1=C1∑j=1N∑k=1Kxijk,where C1 is the fixed use cost of vehicles; if vehicle k is driven directly from the customer i to customer j, then xijk=1; otherwise, it is 0; N is the number of customers served by the logistics distribution center and the number of refrigerated trucks delivered by K.


#### 3.2.2. Transportation Cost of Vehicles

The transportation cost of vehicles mainly refers to the cost of fuel consumption, which is usually proportional to the mileage of vehicles. Therefore, the transportation costΩ2 of vehicles can be expressed as follows:(2)Ω2=C2∑k=1K∑i=0N∑j=0Ndijxijk,where dij is the distance that the vehicle travels directly from customer i to customer j.


#### 3.2.3. Cargo Damage Cost

Unlike ordinary routing optimization, cold chain transportation products are perishable goods. The general distribution routing optimization model considers the damage of goods during the loading and offloading process, which may be the loss caused by collision. However, the cost of damage to fresh agricultural products is mainly caused by temperature changes during cargo transportation and loading and offloading. The refrigerated items are perishable, and factors such as temperature, humidity, oxygen concentration, and product moisture content in the storage environment will affect the changes in product quality. As time and temperature change, the quality of perishable goods will gradually decrease or even lose value. When the quality of the product drops to a certain level, there will be lost costs. The quantity of remaining goods in refrigerated trucks is not only closely related to the needs of customers but also directly related to the loss of goods in the cold chain logistics and distribution process.Variable function of the quality of refrigerated goods:Gt=G0e−∂t, where Gt represents the quality of product at time t and t is the transit time of product; G0 represents the quality of the product when it departs from the distribution center; ∂ represents the corruption rate of the product, and its value is related to the freshness. The characteristics of agricultural products are related to temperature. Because in this study, it is assumed that the ambient temperature of the product does not change during transportation. Therefore, at a certain constant temperature, the corruption rate of fresh agricultural products is constant, and the quality of fresh agricultural products exhibits exponential changes over time [22].Thus, when the delivery vehicle departs from delivery center to customeri, the cargo damage cost Ω31 caused by not opening the door during transportation is as follows:(3)Ω31=∑k=1K∑i=0Nsikpqi1−e∂1tik−t0k,where p is the unit value of fresh agricultural products; qi is the demand of customer i; tik is the time when the refrigerated vehicle k arrives at customer i; t0k is the time when vehicle k departs from distribution center. For the corruption rate of product at a specific temperature: sik is a 0-1 variable; when the vehicle k serves the customer i, then sik=1; otherwise, sik=0.Assuming that the corruption rate at this time is ∂2 (∂2 > ∂1), when the refrigerated truck arrives at customer i, the goods damage cost Ω32 caused by the corruption of fresh agricultural products caused by opening the door during unloading is as follows:(4)Ω32=∑k=1K∑i=0NsikpQin1−e−∂2tsi,where Qin is the remaining product weight when the refrigerated truck leaves customer i. tsi is the service time required by customer i. Therefore, the cargo damage cost Ω3 during the entire distribution process is as follows:(5)Ω3=∑k=1K∑i=0Nsikpqi1−e−∂1tik−t0k+Qin1−e−∂2tsi.


#### 3.2.4. Cooling Cost

The refrigeration cost shall also be considered during transportation. The calculation of costΩ41 includes the energy consumed by the vehicle to maintain low temperature during transportation and the cost of additional energy provided by the refrigeration system during unloading. The calculation is as follows:(6)Ω41=Ce∑k=1K∑i=0N∑j=0Nxijktijk,where Ce is the cost of using refrigerated trucks in the transportation process.During the offloading process, cooling costΩ42 of refrigerated trucks can be expressed as follows:(7)Ω42=Ce′∑k=1K∑j=0Nsjkωj,where Ce′ is the cost of using the refrigerated truck when offloading; ωj is the additional loss cost of customer j.Therefore, the cooling cost of the whole process is as follows:(8)Ω4=Ce∑k=1K∑i=0N∑j=0Nxijktijk+Ce′∑k=1K∑j=0Nsjkωj.


#### 3.2.5. Carbon Emissions Cost

Carbon emissions = fuel consumption × CO2 emission coefficient. Among them, fuel consumption is related to transportation distance and vehicle load capacity. The fuel consumption per unit distance ϑ can be expressed as a linear function depending on the truck load X. If the total vehicle weight is divided into vehicle weight Q0 and load X, the fuel consumption per unit distance ϑX is as follows:(9)ϑX=aQ0+X+b,where a and b are the fuel consumption coefficients.Suppose the maximum cargo capacity of vehicles isQmax, the fuel consumption per unit distance when fully loaded is ϑ∗, and the fuel consumption per unit distance when empty is ϑ0, which is expressed as follows:(10)ϑ0=aQ0+b,ϑ∗=aQ0+Qmax+b,a=ϑ∗−ϑ0Qmax.Thus, the fuel consumption per unit distanceϑX can be expressed as follows:(11)ϑX=ϑ0+ϑ∗−ϑ0QmaxX.Therefore, in the distribution process of fresh agricultural products, if the goods ofQij are transported from customer i to customer j, the carbon emissions generated when driving between i,j can be expressed as follows:(12)E1=τϑQijdij,where τ is CO2 emission coefficient; Qij represents the load capacity of vehicles when the vehicle is directly driven from customer i to customer j; ϑQij represents the fuel consumption per unit distance when the vehicle is traveling between i,j and the load capacity is Qij. The carbon emission cost = carbon tax × carbon emission amount. If the carbon tax is set as C3, then the total carbon emission cost Ω5 in distribution process is as follows:(13)Ω5=C3τ∑k=1K∑i=0N∑j=0NdijxijkϑQij.


#### 3.2.6. Penalty Cost

In the process of distribution, if the goods cannot be delivered within the required time, a certain penalty cost will be incurredΩ6:(14)Ω6=β1∑i=1NmaxLi′−ti,0+β2∑i=1Nmaxti−Li″,0,where β1 and β2 are the penalty cost per unit time when the vehicle arrives at customer i before time Li′ and time Li″; Li′,Li″ is the service time window required by customer i.In summary, the distribution routing optimization model for fresh agricultural products considering carbon emissions is as follows:(15)minΩ=C1∑j=1N∑k=1Kx0jk+C2∑k=1K∑i=0N∑j=0Ndijxijk+∑k=1K∑i=0Nsikpqi1−e−∂1tik−t0k+Qin1−e−∂2tsi+C3τ∑k=1K∑i=0N∑j=0NdijxijkϑQij+Ce∑k=1K∑i=0N∑j=0Nxijktijk+Ce′∑k=1K∑j=0Nsjkωjβ1∑i=1NmaxLi′−ti,0+β2∑i=1Nmaxti−Li″,0.The constraints are as follows:①
Constraints on vehicle carrying capacity:∑i=1Nqisik≤Qk,∀k.②
Ensure that each customer is served by only one vehicle:∑k=1Ksik=1,∀i.③
For any customer, only one vehicle can arrive and depart once:(16)∑i=0Nxijk=sjk,∀j,k,∑j=0Nxijk=sik,∀i,k.④
Elimination condition of the secondary loop:(17)∑i,j∈S×Sxijk≤S−1,S⊆1,2,…,N.⑤
Continuity of the distribution process:tj=ti+tij.⑥
Variable constraints:xijk=0or1,∀i,j,kyik=0or1,∀i,k.


## 3.1. Problem Description

With the expansion scale of cities and the relative dispersion of commercial centers, there is an increasing demand for urban multiregional distribution and chain supermarket distribution. How to realize the effective integration of logistics resources and reduce the cost of urban logistics distribution has become an urgent problem to be solved in current urban logistics distribution [18]. The urban multiregional distribution problem is a typical MDVRP problem involving multiple distribution centers participating in the fulfillment of distribution requirements. The main difference between it and the routing optimization problem of a single distribution center is how different distribution centers in MDVRP jointly participate in the distribution process [19]. At present, large-scale vehicles in cities are affected by traffic restriction policy, and most of the enterprises’ logistics distribution adopts hierarchical distribution. Set up large warehouses on the outskirts of the city and deliver goods to various distribution centers through unrestricted small and medium-sized vehicles. The distribution centers adopt a partitioned distribution model to be responsible for customer service in each area, as shown in Figure 1.Figure 1
Schematic diagram of VRP in multiple distribution centers.MDVRP is generally described as follows: in the logistics network, multiple distribution centers are responsible for the distribution services of customers, and the locations of distribution centers and customers in the static MDVRP are known. At the same time, it is assumed that the goods in the warehouse and the distribution center can meet the needs of all customers, and the goods of different customer needs can be mixed with each other. It is required that the demand of each customer does not exceed the maximum load capacity of vehicles and that only one vehicle is responsible, and the demand is not allowed to be split [20, 21]. It is assumed that the maximum load of vehicles is known and overloading is not allowed. The maximum travel distance for one delivery is known, and violations are not allowed. The distribution distance between the distribution center and customers and between the customers is known. It is required that the vehicles start from the distribution center when the distribution task starts, and it is required to return to the original distribution center after the task is completed.


## 3.2. Model Establishment

### 3.2.1. Fixed Cost of Vehicles

The fixed cost of using a vehicle is usually constant, including the fixed loss of vehicles and the cost related to the use of vehicles, such as the driver’s salary. It has nothing to do with the mileage of vehicles and the number of customers served. The fixed costΩ1 of vehicles is calculated as follows:(1)Ω1=C1∑j=1N∑k=1Kxijk,where C1 is the fixed use cost of vehicles; if vehicle k is driven directly from the customer i to customer j, then xijk=1; otherwise, it is 0; N is the number of customers served by the logistics distribution center and the number of refrigerated trucks delivered by K.


### 3.2.2. Transportation Cost of Vehicles

The transportation cost of vehicles mainly refers to the cost of fuel consumption, which is usually proportional to the mileage of vehicles. Therefore, the transportation costΩ2 of vehicles can be expressed as follows:(2)Ω2=C2∑k=1K∑i=0N∑j=0Ndijxijk,where dij is the distance that the vehicle travels directly from customer i to customer j.


### 3.2.3. Cargo Damage Cost

Unlike ordinary routing optimization, cold chain transportation products are perishable goods. The general distribution routing optimization model considers the damage of goods during the loading and offloading process, which may be the loss caused by collision. However, the cost of damage to fresh agricultural products is mainly caused by temperature changes during cargo transportation and loading and offloading. The refrigerated items are perishable, and factors such as temperature, humidity, oxygen concentration, and product moisture content in the storage environment will affect the changes in product quality. As time and temperature change, the quality of perishable goods will gradually decrease or even lose value. When the quality of the product drops to a certain level, there will be lost costs. The quantity of remaining goods in refrigerated trucks is not only closely related to the needs of customers but also directly related to the loss of goods in the cold chain logistics and distribution process.Variable function of the quality of refrigerated goods:Gt=G0e−∂t, where Gt represents the quality of product at time t and t is the transit time of product; G0 represents the quality of the product when it departs from the distribution center; ∂ represents the corruption rate of the product, and its value is related to the freshness. The characteristics of agricultural products are related to temperature. Because in this study, it is assumed that the ambient temperature of the product does not change during transportation. Therefore, at a certain constant temperature, the corruption rate of fresh agricultural products is constant, and the quality of fresh agricultural products exhibits exponential changes over time [22].Thus, when the delivery vehicle departs from delivery center to customeri, the cargo damage cost Ω31 caused by not opening the door during transportation is as follows:(3)Ω31=∑k=1K∑i=0Nsikpqi1−e∂1tik−t0k,where p is the unit value of fresh agricultural products; qi is the demand of customer i; tik is the time when the refrigerated vehicle k arrives at customer i; t0k is the time when vehicle k departs from distribution center. For the corruption rate of product at a specific temperature: sik is a 0-1 variable; when the vehicle k serves the customer i, then sik=1; otherwise, sik=0.Assuming that the corruption rate at this time is ∂2 (∂2 > ∂1), when the refrigerated truck arrives at customer i, the goods damage cost Ω32 caused by the corruption of fresh agricultural products caused by opening the door during unloading is as follows:(4)Ω32=∑k=1K∑i=0NsikpQin1−e−∂2tsi,where Qin is the remaining product weight when the refrigerated truck leaves customer i. tsi is the service time required by customer i. Therefore, the cargo damage cost Ω3 during the entire distribution process is as follows:(5)Ω3=∑k=1K∑i=0Nsikpqi1−e−∂1tik−t0k+Qin1−e−∂2tsi.


### 3.2.4. Cooling Cost

The refrigeration cost shall also be considered during transportation. The calculation of costΩ41 includes the energy consumed by the vehicle to maintain low temperature during transportation and the cost of additional energy provided by the refrigeration system during unloading. The calculation is as follows:(6)Ω41=Ce∑k=1K∑i=0N∑j=0Nxijktijk,where Ce is the cost of using refrigerated trucks in the transportation process.During the offloading process, cooling costΩ42 of refrigerated trucks can be expressed as follows:(7)Ω42=Ce′∑k=1K∑j=0Nsjkωj,where Ce′ is the cost of using the refrigerated truck when offloading; ωj is the additional loss cost of customer j.Therefore, the cooling cost of the whole process is as follows:(8)Ω4=Ce∑k=1K∑i=0N∑j=0Nxijktijk+Ce′∑k=1K∑j=0Nsjkωj.


### 3.2.5. Carbon Emissions Cost

Carbon emissions = fuel consumption × CO2 emission coefficient. Among them, fuel consumption is related to transportation distance and vehicle load capacity. The fuel consumption per unit distance ϑ can be expressed as a linear function depending on the truck load X. If the total vehicle weight is divided into vehicle weight Q0 and load X, the fuel consumption per unit distance ϑX is as follows:(9)ϑX=aQ0+X+b,where a and b are the fuel consumption coefficients.Suppose the maximum cargo capacity of vehicles isQmax, the fuel consumption per unit distance when fully loaded is ϑ∗, and the fuel consumption per unit distance when empty is ϑ0, which is expressed as follows:(10)ϑ0=aQ0+b,ϑ∗=aQ0+Qmax+b,a=ϑ∗−ϑ0Qmax.Thus, the fuel consumption per unit distanceϑX can be expressed as follows:(11)ϑX=ϑ0+ϑ∗−ϑ0QmaxX.Therefore, in the distribution process of fresh agricultural products, if the goods ofQij are transported from customer i to customer j, the carbon emissions generated when driving between i,j can be expressed as follows:(12)E1=τϑQijdij,where τ is CO2 emission coefficient; Qij represents the load capacity of vehicles when the vehicle is directly driven from customer i to customer j; ϑQij represents the fuel consumption per unit distance when the vehicle is traveling between i,j and the load capacity is Qij. The carbon emission cost = carbon tax × carbon emission amount. If the carbon tax is set as C3, then the total carbon emission cost Ω5 in distribution process is as follows:(13)Ω5=C3τ∑k=1K∑i=0N∑j=0NdijxijkϑQij.


### 3.2.6. Penalty Cost

In the process of distribution, if the goods cannot be delivered within the required time, a certain penalty cost will be incurredΩ6:(14)Ω6=β1∑i=1NmaxLi′−ti,0+β2∑i=1Nmaxti−Li″,0,where β1 and β2 are the penalty cost per unit time when the vehicle arrives at customer i before time Li′ and time Li″; Li′,Li″ is the service time window required by customer i.In summary, the distribution routing optimization model for fresh agricultural products considering carbon emissions is as follows:(15)minΩ=C1∑j=1N∑k=1Kx0jk+C2∑k=1K∑i=0N∑j=0Ndijxijk+∑k=1K∑i=0Nsikpqi1−e−∂1tik−t0k+Qin1−e−∂2tsi+C3τ∑k=1K∑i=0N∑j=0NdijxijkϑQij+Ce∑k=1K∑i=0N∑j=0Nxijktijk+Ce′∑k=1K∑j=0Nsjkωjβ1∑i=1NmaxLi′−ti,0+β2∑i=1Nmaxti−Li″,0.The constraints are as follows:①
Constraints on vehicle carrying capacity:∑i=1Nqisik≤Qk,∀k.②
Ensure that each customer is served by only one vehicle:∑k=1Ksik=1,∀i.③
For any customer, only one vehicle can arrive and depart once:(16)∑i=0Nxijk=sjk,∀j,k,∑j=0Nxijk=sik,∀i,k.④
Elimination condition of the secondary loop:(17)∑i,j∈S×Sxijk≤S−1,S⊆1,2,…,N.⑤
Continuity of the distribution process:tj=ti+tij.⑥
Variable constraints:xijk=0or1,∀i,j,kyik=0or1,∀i,k.


## 3.2.1. Fixed Cost of Vehicles

The fixed cost of using a vehicle is usually constant, including the fixed loss of vehicles and the cost related to the use of vehicles, such as the driver’s salary. It has nothing to do with the mileage of vehicles and the number of customers served. The fixed costΩ1 of vehicles is calculated as follows:(1)Ω1=C1∑j=1N∑k=1Kxijk,where C1 is the fixed use cost of vehicles; if vehicle k is driven directly from the customer i to customer j, then xijk=1; otherwise, it is 0; N is the number of customers served by the logistics distribution center and the number of refrigerated trucks delivered by K.


## 3.2.2. Transportation Cost of Vehicles

The transportation cost of vehicles mainly refers to the cost of fuel consumption, which is usually proportional to the mileage of vehicles. Therefore, the transportation costΩ2 of vehicles can be expressed as follows:(2)Ω2=C2∑k=1K∑i=0N∑j=0Ndijxijk,where dij is the distance that the vehicle travels directly from customer i to customer j.


## 3.2.3. Cargo Damage Cost

Unlike ordinary routing optimization, cold chain transportation products are perishable goods. The general distribution routing optimization model considers the damage of goods during the loading and offloading process, which may be the loss caused by collision. However, the cost of damage to fresh agricultural products is mainly caused by temperature changes during cargo transportation and loading and offloading. The refrigerated items are perishable, and factors such as temperature, humidity, oxygen concentration, and product moisture content in the storage environment will affect the changes in product quality. As time and temperature change, the quality of perishable goods will gradually decrease or even lose value. When the quality of the product drops to a certain level, there will be lost costs. The quantity of remaining goods in refrigerated trucks is not only closely related to the needs of customers but also directly related to the loss of goods in the cold chain logistics and distribution process.Variable function of the quality of refrigerated goods:Gt=G0e−∂t, where Gt represents the quality of product at time t and t is the transit time of product; G0 represents the quality of the product when it departs from the distribution center; ∂ represents the corruption rate of the product, and its value is related to the freshness. The characteristics of agricultural products are related to temperature. Because in this study, it is assumed that the ambient temperature of the product does not change during transportation. Therefore, at a certain constant temperature, the corruption rate of fresh agricultural products is constant, and the quality of fresh agricultural products exhibits exponential changes over time [22].Thus, when the delivery vehicle departs from delivery center to customeri, the cargo damage cost Ω31 caused by not opening the door during transportation is as follows:(3)Ω31=∑k=1K∑i=0Nsikpqi1−e∂1tik−t0k,where p is the unit value of fresh agricultural products; qi is the demand of customer i; tik is the time when the refrigerated vehicle k arrives at customer i; t0k is the time when vehicle k departs from distribution center. For the corruption rate of product at a specific temperature: sik is a 0-1 variable; when the vehicle k serves the customer i, then sik=1; otherwise, sik=0.Assuming that the corruption rate at this time is ∂2 (∂2 > ∂1), when the refrigerated truck arrives at customer i, the goods damage cost Ω32 caused by the corruption of fresh agricultural products caused by opening the door during unloading is as follows:(4)Ω32=∑k=1K∑i=0NsikpQin1−e−∂2tsi,where Qin is the remaining product weight when the refrigerated truck leaves customer i. tsi is the service time required by customer i. Therefore, the cargo damage cost Ω3 during the entire distribution process is as follows:(5)Ω3=∑k=1K∑i=0Nsikpqi1−e−∂1tik−t0k+Qin1−e−∂2tsi.


## 3.2.4. Cooling Cost

The refrigeration cost shall also be considered during transportation. The calculation of costΩ41 includes the energy consumed by the vehicle to maintain low temperature during transportation and the cost of additional energy provided by the refrigeration system during unloading. The calculation is as follows:(6)Ω41=Ce∑k=1K∑i=0N∑j=0Nxijktijk,where Ce is the cost of using refrigerated trucks in the transportation process.During the offloading process, cooling costΩ42 of refrigerated trucks can be expressed as follows:(7)Ω42=Ce′∑k=1K∑j=0Nsjkωj,where Ce′ is the cost of using the refrigerated truck when offloading; ωj is the additional loss cost of customer j.Therefore, the cooling cost of the whole process is as follows:(8)Ω4=Ce∑k=1K∑i=0N∑j=0Nxijktijk+Ce′∑k=1K∑j=0Nsjkωj.


## 3.2.5. Carbon Emissions Cost

Carbon emissions = fuel consumption × CO2 emission coefficient. Among them, fuel consumption is related to transportation distance and vehicle load capacity. The fuel consumption per unit distance ϑ can be expressed as a linear function depending on the truck load X. If the total vehicle weight is divided into vehicle weight Q0 and load X, the fuel consumption per unit distance ϑX is as follows:(9)ϑX=aQ0+X+b,where a and b are the fuel consumption coefficients.Suppose the maximum cargo capacity of vehicles isQmax, the fuel consumption per unit distance when fully loaded is ϑ∗, and the fuel consumption per unit distance when empty is ϑ0, which is expressed as follows:(10)ϑ0=aQ0+b,ϑ∗=aQ0+Qmax+b,a=ϑ∗−ϑ0Qmax.Thus, the fuel consumption per unit distanceϑX can be expressed as follows:(11)ϑX=ϑ0+ϑ∗−ϑ0QmaxX.Therefore, in the distribution process of fresh agricultural products, if the goods ofQij are transported from customer i to customer j, the carbon emissions generated when driving between i,j can be expressed as follows:(12)E1=τϑQijdij,where τ is CO2 emission coefficient; Qij represents the load capacity of vehicles when the vehicle is directly driven from customer i to customer j; ϑQij represents the fuel consumption per unit distance when the vehicle is traveling between i,j and the load capacity is Qij. The carbon emission cost = carbon tax × carbon emission amount. If the carbon tax is set as C3, then the total carbon emission cost Ω5 in distribution process is as follows:(13)Ω5=C3τ∑k=1K∑i=0N∑j=0NdijxijkϑQij.


## 3.2.6. Penalty Cost

In the process of distribution, if the goods cannot be delivered within the required time, a certain penalty cost will be incurredΩ6:(14)Ω6=β1∑i=1NmaxLi′−ti,0+β2∑i=1Nmaxti−Li″,0,where β1 and β2 are the penalty cost per unit time when the vehicle arrives at customer i before time Li′ and time Li″; Li′,Li″ is the service time window required by customer i.In summary, the distribution routing optimization model for fresh agricultural products considering carbon emissions is as follows:(15)minΩ=C1∑j=1N∑k=1Kx0jk+C2∑k=1K∑i=0N∑j=0Ndijxijk+∑k=1K∑i=0Nsikpqi1−e−∂1tik−t0k+Qin1−e−∂2tsi+C3τ∑k=1K∑i=0N∑j=0NdijxijkϑQij+Ce∑k=1K∑i=0N∑j=0Nxijktijk+Ce′∑k=1K∑j=0Nsjkωjβ1∑i=1NmaxLi′−ti,0+β2∑i=1Nmaxti−Li″,0.The constraints are as follows:①
Constraints on vehicle carrying capacity:∑i=1Nqisik≤Qk,∀k.②
Ensure that each customer is served by only one vehicle:∑k=1Ksik=1,∀i.③
For any customer, only one vehicle can arrive and depart once:(16)∑i=0Nxijk=sjk,∀j,k,∑j=0Nxijk=sik,∀i,k.④
Elimination condition of the secondary loop:(17)∑i,j∈S×Sxijk≤S−1,S⊆1,2,…,N.⑤
Continuity of the distribution process:tj=ti+tij.⑥
Variable constraints:xijk=0or1,∀i,j,kyik=0or1,∀i,k.


## 4. Algorithm Design

### 4.1. Genetic Algorithm

Genetic algorithm (GA) can simulate the principle of biological evolution in nature and draw inspiration from it. It is based on the basic idea of a natural gene screening mechanism to form a set of simulated evolutionary algorithms to solve complex problems in reality [23]. GA is widely used in many disciplines and has a core position in intelligent algorithms.GA includes five elements: parameter coding, initial population setting, fitness function setting, genetic operation setting, and control parameter setting. Its structure is shown in Figure2.Figure 2
Population structure evolution.GA is different from precise algorithms and enumeration methods. It is a global random search algorithm, which often uses genetic principles and biological evolution theories as references. GA converts the digital information that needs to be processed into chromosomes in the form of encoding and then simulates the principle of species evolution. Exchange the information on the chromosomes to find the optimal solution [24, 25]. At the same time, GA has the characteristics of randomness, which aims to prevent local premature convergence. According to chromosomal fitness, the highly adaptable individuals are selected to reconstitute a new population and enter the inheritance of the next generation.GA is widely used in practical problems. As a general algorithm, its advantages include wide application fields; complex problems in real life can be solved by GA. It has scalability and is compatible with various problems in various fields, which is conducive to solving optimization problems. However, there are also certain shortcomings: the operation efficiency is generally lower than other algorithms, the coding does not have a unified specification, and the representation is not accurate enough. The algorithm lacks effective quantitative analysis, especially in terms of accuracy and reliability.


### 4.2. Improved GA for Distribution Routing Optimization

#### 4.2.1. Encoding

Instead of using binary encoding, natural number encoding is used to encode chromosomes [26]. The main decision variables in this study are the number of vehicles dispatched and the number of customers served by each vehicle and the order of service. Therefore, the chromosomes are coded using the corresponding arrangement of vehicles and customers, and each chromosome represents a solution. The specific idea is to use N1∼N natural numbers without repetition to represent the number of customers, and the arrangement of the natural numbers represents the order of serving customers. Then use any repeatable natural number arrangement among N1∼R to represent the vehicle arrangement. Corresponding between the two is a solution, which corresponds to a distribution plan.For example, use natural numbers 1 to 8 to represent 8 customers, and then randomly generate a full array of 8 customers (2 5 4 6 7 3 1 8). Then generate the corresponding vehicle arrangement (1 1 1 2 2 3 3 3) based on the vehicle load and other constraints, which means that the first car serves customers 2, 5, and 4, the second car serves customers 6 and 7, and the third vehicle serves customers 3, 1, and 8. The specific method for generating the corresponding vehicle arrangement is as follows: First, the full arrangement ofN customers is expressed as L=L1,L2,…,Ln, and customer L1∼Lr1 is assigned to vehicle 1, r1∈r1|∑i=1r1qi≤Q,∑i=1r1+1qi≥Q,1≤r1≤S, where qi is the demand of customer i and Q is the maximum load of vehicles. Then, the customer Lr1+1∼Lr2 is handed over to the transportation vehicle 2 for service and so on to assign tasks.


#### 4.2.2. Population Initialization

Chaos research first originated in the fields of physics and mathematics, and then under the background of the development of computer technology, its research scope was expanded to the fields of natural science and social science. Since the chaos optimization algorithm has good search traversal, the chaotic state can be introduced into the optimization variable to expand the traversal range so that the chaotic variable can be quickly and conveniently searched. A large number of reference studies have shown that the Chaotic GA (CGA) has a very good complementary performance in solving large-scale uncertain selection problems, is not easy to fall into a local optimum, increases the diversity of the population, and makes the algorithm closer to real organisms. Evolutionary process. Combining the ergodicity of chaotic search and the reproductive characteristics of GA, chaotic search is introduced into optimization variables. Encode chaotic variables into chromosomes to form individuals to participate in the genetic process and finally get the best individual.The chaos technology is introduced into the population initialization; the purpose is to make the quality of the initial population of GA better and the calculation efficiency higher. The selection of the initial population is very important and directly affects the entire genetic process. The introduction of chaotic initial conditions will change the current situation of poor local optimization ability and easy “premature maturity.” Commonly used chaos models mainly include logistic mapping and Tent mapping.Logistic mapping:(18)xn=ϖxn1−xn;x∈0,1.Tent mapping:(19)xn+1=2xn,0≤xn≤0.521−xn0.5≤xn≤1n=0,1,2,….The difference between the two lies in the probability density distribution, and the logistic mapping density function has an uneven distribution. And the mapping is in the case of 30,000 iterations in the [0, 1] interval, and the value density is higher in the [0, 0.1] and [0.9, 1] intervals, reaching 1,900 times. In the interval of [0.1, 0.9], the average value is only 200 times, and the uniform distribution is poor. Similarly, when the Tent map is iterated 30000 times in the [0, 1] interval, the highest value is 439 times, the lowest value is 248 times, and the average is 300 times. The distribution ratio is uniform, which can speed up the search speed. Based on this, the Tent map is selected to generate the turbidity sequence.In the initialization process of Tent chaos,U random integers with a value range of [0, n] are first generated and combined into a chromosome set a1,a2,…,ai,…,aU. Then randomly generate U random numbers with a value range of [0,1], use the Tent mapping expression to chaotically disturb them, generate chaotic variables, and form a chaotic sequence Xij. Then map the chaotic variable to the value range of the optimized variable according to the following formula and round it:(20)gij=1+n−1Xij,i=1,2,…,U;j=1,2,…,k.Obtain the chaotic chromosomeYi=gi1,gi2,…,giki, delete the same genes in them, delete the adjacent 0 when multiple 0s are adjacent, and keep the adjacent one. Calculate the fitness of each chromosome Yi, and select the top TOP M individuals as the initial population. The initial population size in this article is 100.


#### 4.2.3. Fitness Function

Fitness function is not only an important basis for genetic operation but also a standard to measure the quality of individuals. The larger the value of the fitness function, the better. However, in practical problems, the objective function is not always nonnegative, sometimes the maximum value and sometimes the minimum value.The greater the fitness of an individual, the better the performance of the individual. However, the objective function of the green and low-carbon cold chain distribution routing optimization problem is the minimum cost. Therefore, the objective function needs to be converted into fitness, and the reciprocal of the objective function can be expressed as the fitness function. And because the optimization problem under the time window constraint is established, the time constraint can be dealt with by adding a penalty function. If the solution corresponding to an individual violates the time constraint, a certain penalty will be given to make it have a smaller fitness, and the probability of being inherited will also become smaller. Therefore, the constructed fitness function is expressed as the inverse of the sum of the objective cost function and the penalty cost of violating the time window, expressed as follows:(21)fi=1zi+ψi,where fi represents the fitness corresponding to individual i, zi represents the objective function value corresponding to individual i, and ψi represents the penalty cost of individual i for violating the time window.


#### 4.2.4. Selection Operator

The selection operation determines the probability of being selected according to the size of the fitness value. It is based on the evaluation of the individual’s fitness. The higher the evaluation is, the easier it is to be selected and inherited. The probability of individual selection is calculated as follows:(22)Pk=fk∑i=1nfi.First, calculate the fitness of each body in each generation; sum and then calculate the total value of the individual and the fitness, that is, the probability of the individual’s selection being inherited. Sort all individual selection probabilities. The higher the probability, the better the individual performance, and the higher the ranking, the easier it is to be copied to the next generation.After the selection operation, sort according to the fitness value, and select the firstm−1 chromosomes to perform turbidity perturbation again according to equation (20). Calculate the fitness of each perturbed chromosome, inherit the chromosome individual with the largest fitness value directly to the next generation, and perform the next crossover and mutation operations on the remaining chromosomes. This makes the outstanding individuals after the selection operation chaotic again, increases the diversity of the population, avoids premature maturity, and accelerates the convergence speed of the algorithm.


#### 4.2.5. Crossover Operator

The cross-reorganization is carried out using the cyclic crossover rule, and the specific operation process is shown in Figure3.Figure 3
Cyclic crossover operator.Firstly, find a cycle based on the corresponding customer positions in the two individuals of the parent. Then, copy the customer in the cycle position of an individual in the parent to a descendant. Finally, delete the customer in the recurring position of another individual in the parent, and copy the customer in the remaining nonrecurring position directly to the remaining position of the offspring.


#### 4.2.6. Mutation Operator

The mutation operation is to change the value of some genes of an individual. The use of mutation operators in GA can improve the algorithm’s ability to solve problems, speed up the generation of optimal solutions, and make the population more diverse and reduce the probability of premature maturity. According to the characteristics of the optimized problem, the use of inversion mutation is also called inverse transformation. Examples are as follows:(1)
Randomly select 2 points from the individual. For example,temp = 589437612; randomly select the 3rd and 6th points; that is, temp=589437612; “∥” represents the crossover sequence;(2)
Perform the reverse transformation operation to obtain the mutated individualtemp1=587349621. The mutation operation is shown in Figure 4.Figure 4
Mutation operation.


#### 4.2.7. Termination Criteria

The termination criterion is a condition for judging whether to stop the operation, and the selected termination condition is to reach the preset 2000 generations of evolution. If it reaches the evolutionary algebra, then stop the evolution. In the last-generation population, the distribution path corresponding to the individual with the best fitness value is the optimal solution to the green cold chain logistics distribution routing optimization problem.


## 4.1. Genetic Algorithm

Genetic algorithm (GA) can simulate the principle of biological evolution in nature and draw inspiration from it. It is based on the basic idea of a natural gene screening mechanism to form a set of simulated evolutionary algorithms to solve complex problems in reality [23]. GA is widely used in many disciplines and has a core position in intelligent algorithms.GA includes five elements: parameter coding, initial population setting, fitness function setting, genetic operation setting, and control parameter setting. Its structure is shown in Figure2.Figure 2
Population structure evolution.GA is different from precise algorithms and enumeration methods. It is a global random search algorithm, which often uses genetic principles and biological evolution theories as references. GA converts the digital information that needs to be processed into chromosomes in the form of encoding and then simulates the principle of species evolution. Exchange the information on the chromosomes to find the optimal solution [24, 25]. At the same time, GA has the characteristics of randomness, which aims to prevent local premature convergence. According to chromosomal fitness, the highly adaptable individuals are selected to reconstitute a new population and enter the inheritance of the next generation.GA is widely used in practical problems. As a general algorithm, its advantages include wide application fields; complex problems in real life can be solved by GA. It has scalability and is compatible with various problems in various fields, which is conducive to solving optimization problems. However, there are also certain shortcomings: the operation efficiency is generally lower than other algorithms, the coding does not have a unified specification, and the representation is not accurate enough. The algorithm lacks effective quantitative analysis, especially in terms of accuracy and reliability.


## 4.2. Improved GA for Distribution Routing Optimization

### 4.2.1. Encoding

Instead of using binary encoding, natural number encoding is used to encode chromosomes [26]. The main decision variables in this study are the number of vehicles dispatched and the number of customers served by each vehicle and the order of service. Therefore, the chromosomes are coded using the corresponding arrangement of vehicles and customers, and each chromosome represents a solution. The specific idea is to use N1∼N natural numbers without repetition to represent the number of customers, and the arrangement of the natural numbers represents the order of serving customers. Then use any repeatable natural number arrangement among N1∼R to represent the vehicle arrangement. Corresponding between the two is a solution, which corresponds to a distribution plan.For example, use natural numbers 1 to 8 to represent 8 customers, and then randomly generate a full array of 8 customers (2 5 4 6 7 3 1 8). Then generate the corresponding vehicle arrangement (1 1 1 2 2 3 3 3) based on the vehicle load and other constraints, which means that the first car serves customers 2, 5, and 4, the second car serves customers 6 and 7, and the third vehicle serves customers 3, 1, and 8. The specific method for generating the corresponding vehicle arrangement is as follows: First, the full arrangement ofN customers is expressed as L=L1,L2,…,Ln, and customer L1∼Lr1 is assigned to vehicle 1, r1∈r1|∑i=1r1qi≤Q,∑i=1r1+1qi≥Q,1≤r1≤S, where qi is the demand of customer i and Q is the maximum load of vehicles. Then, the customer Lr1+1∼Lr2 is handed over to the transportation vehicle 2 for service and so on to assign tasks.


### 4.2.2. Population Initialization

Chaos research first originated in the fields of physics and mathematics, and then under the background of the development of computer technology, its research scope was expanded to the fields of natural science and social science. Since the chaos optimization algorithm has good search traversal, the chaotic state can be introduced into the optimization variable to expand the traversal range so that the chaotic variable can be quickly and conveniently searched. A large number of reference studies have shown that the Chaotic GA (CGA) has a very good complementary performance in solving large-scale uncertain selection problems, is not easy to fall into a local optimum, increases the diversity of the population, and makes the algorithm closer to real organisms. Evolutionary process. Combining the ergodicity of chaotic search and the reproductive characteristics of GA, chaotic search is introduced into optimization variables. Encode chaotic variables into chromosomes to form individuals to participate in the genetic process and finally get the best individual.The chaos technology is introduced into the population initialization; the purpose is to make the quality of the initial population of GA better and the calculation efficiency higher. The selection of the initial population is very important and directly affects the entire genetic process. The introduction of chaotic initial conditions will change the current situation of poor local optimization ability and easy “premature maturity.” Commonly used chaos models mainly include logistic mapping and Tent mapping.Logistic mapping:(18)xn=ϖxn1−xn;x∈0,1.Tent mapping:(19)xn+1=2xn,0≤xn≤0.521−xn0.5≤xn≤1n=0,1,2,….The difference between the two lies in the probability density distribution, and the logistic mapping density function has an uneven distribution. And the mapping is in the case of 30,000 iterations in the [0, 1] interval, and the value density is higher in the [0, 0.1] and [0.9, 1] intervals, reaching 1,900 times. In the interval of [0.1, 0.9], the average value is only 200 times, and the uniform distribution is poor. Similarly, when the Tent map is iterated 30000 times in the [0, 1] interval, the highest value is 439 times, the lowest value is 248 times, and the average is 300 times. The distribution ratio is uniform, which can speed up the search speed. Based on this, the Tent map is selected to generate the turbidity sequence.In the initialization process of Tent chaos,U random integers with a value range of [0, n] are first generated and combined into a chromosome set a1,a2,…,ai,…,aU. Then randomly generate U random numbers with a value range of [0,1], use the Tent mapping expression to chaotically disturb them, generate chaotic variables, and form a chaotic sequence Xij. Then map the chaotic variable to the value range of the optimized variable according to the following formula and round it:(20)gij=1+n−1Xij,i=1,2,…,U;j=1,2,…,k.Obtain the chaotic chromosomeYi=gi1,gi2,…,giki, delete the same genes in them, delete the adjacent 0 when multiple 0s are adjacent, and keep the adjacent one. Calculate the fitness of each chromosome Yi, and select the top TOP M individuals as the initial population. The initial population size in this article is 100.


### 4.2.3. Fitness Function

Fitness function is not only an important basis for genetic operation but also a standard to measure the quality of individuals. The larger the value of the fitness function, the better. However, in practical problems, the objective function is not always nonnegative, sometimes the maximum value and sometimes the minimum value.The greater the fitness of an individual, the better the performance of the individual. However, the objective function of the green and low-carbon cold chain distribution routing optimization problem is the minimum cost. Therefore, the objective function needs to be converted into fitness, and the reciprocal of the objective function can be expressed as the fitness function. And because the optimization problem under the time window constraint is established, the time constraint can be dealt with by adding a penalty function. If the solution corresponding to an individual violates the time constraint, a certain penalty will be given to make it have a smaller fitness, and the probability of being inherited will also become smaller. Therefore, the constructed fitness function is expressed as the inverse of the sum of the objective cost function and the penalty cost of violating the time window, expressed as follows:(21)fi=1zi+ψi,where fi represents the fitness corresponding to individual i, zi represents the objective function value corresponding to individual i, and ψi represents the penalty cost of individual i for violating the time window.


### 4.2.4. Selection Operator

The selection operation determines the probability of being selected according to the size of the fitness value. It is based on the evaluation of the individual’s fitness. The higher the evaluation is, the easier it is to be selected and inherited. The probability of individual selection is calculated as follows:(22)Pk=fk∑i=1nfi.First, calculate the fitness of each body in each generation; sum and then calculate the total value of the individual and the fitness, that is, the probability of the individual’s selection being inherited. Sort all individual selection probabilities. The higher the probability, the better the individual performance, and the higher the ranking, the easier it is to be copied to the next generation.After the selection operation, sort according to the fitness value, and select the firstm−1 chromosomes to perform turbidity perturbation again according to equation (20). Calculate the fitness of each perturbed chromosome, inherit the chromosome individual with the largest fitness value directly to the next generation, and perform the next crossover and mutation operations on the remaining chromosomes. This makes the outstanding individuals after the selection operation chaotic again, increases the diversity of the population, avoids premature maturity, and accelerates the convergence speed of the algorithm.


### 4.2.5. Crossover Operator

The cross-reorganization is carried out using the cyclic crossover rule, and the specific operation process is shown in Figure3.Figure 3
Cyclic crossover operator.Firstly, find a cycle based on the corresponding customer positions in the two individuals of the parent. Then, copy the customer in the cycle position of an individual in the parent to a descendant. Finally, delete the customer in the recurring position of another individual in the parent, and copy the customer in the remaining nonrecurring position directly to the remaining position of the offspring.


### 4.2.6. Mutation Operator

The mutation operation is to change the value of some genes of an individual. The use of mutation operators in GA can improve the algorithm’s ability to solve problems, speed up the generation of optimal solutions, and make the population more diverse and reduce the probability of premature maturity. According to the characteristics of the optimized problem, the use of inversion mutation is also called inverse transformation. Examples are as follows:(1)
Randomly select 2 points from the individual. For example,temp = 589437612; randomly select the 3rd and 6th points; that is, temp=589437612; “∥” represents the crossover sequence;(2)
Perform the reverse transformation operation to obtain the mutated individualtemp1=587349621. The mutation operation is shown in Figure 4.Figure 4
Mutation operation.


### 4.2.7. Termination Criteria

The termination criterion is a condition for judging whether to stop the operation, and the selected termination condition is to reach the preset 2000 generations of evolution. If it reaches the evolutionary algebra, then stop the evolution. In the last-generation population, the distribution path corresponding to the individual with the best fitness value is the optimal solution to the green cold chain logistics distribution routing optimization problem.


## 4.2.1. Encoding

Instead of using binary encoding, natural number encoding is used to encode chromosomes [26]. The main decision variables in this study are the number of vehicles dispatched and the number of customers served by each vehicle and the order of service. Therefore, the chromosomes are coded using the corresponding arrangement of vehicles and customers, and each chromosome represents a solution. The specific idea is to use N1∼N natural numbers without repetition to represent the number of customers, and the arrangement of the natural numbers represents the order of serving customers. Then use any repeatable natural number arrangement among N1∼R to represent the vehicle arrangement. Corresponding between the two is a solution, which corresponds to a distribution plan.For example, use natural numbers 1 to 8 to represent 8 customers, and then randomly generate a full array of 8 customers (2 5 4 6 7 3 1 8). Then generate the corresponding vehicle arrangement (1 1 1 2 2 3 3 3) based on the vehicle load and other constraints, which means that the first car serves customers 2, 5, and 4, the second car serves customers 6 and 7, and the third vehicle serves customers 3, 1, and 8. The specific method for generating the corresponding vehicle arrangement is as follows: First, the full arrangement ofN customers is expressed as L=L1,L2,…,Ln, and customer L1∼Lr1 is assigned to vehicle 1, r1∈r1|∑i=1r1qi≤Q,∑i=1r1+1qi≥Q,1≤r1≤S, where qi is the demand of customer i and Q is the maximum load of vehicles. Then, the customer Lr1+1∼Lr2 is handed over to the transportation vehicle 2 for service and so on to assign tasks.


## 4.2.2. Population Initialization

Chaos research first originated in the fields of physics and mathematics, and then under the background of the development of computer technology, its research scope was expanded to the fields of natural science and social science. Since the chaos optimization algorithm has good search traversal, the chaotic state can be introduced into the optimization variable to expand the traversal range so that the chaotic variable can be quickly and conveniently searched. A large number of reference studies have shown that the Chaotic GA (CGA) has a very good complementary performance in solving large-scale uncertain selection problems, is not easy to fall into a local optimum, increases the diversity of the population, and makes the algorithm closer to real organisms. Evolutionary process. Combining the ergodicity of chaotic search and the reproductive characteristics of GA, chaotic search is introduced into optimization variables. Encode chaotic variables into chromosomes to form individuals to participate in the genetic process and finally get the best individual.The chaos technology is introduced into the population initialization; the purpose is to make the quality of the initial population of GA better and the calculation efficiency higher. The selection of the initial population is very important and directly affects the entire genetic process. The introduction of chaotic initial conditions will change the current situation of poor local optimization ability and easy “premature maturity.” Commonly used chaos models mainly include logistic mapping and Tent mapping.Logistic mapping:(18)xn=ϖxn1−xn;x∈0,1.Tent mapping:(19)xn+1=2xn,0≤xn≤0.521−xn0.5≤xn≤1n=0,1,2,….The difference between the two lies in the probability density distribution, and the logistic mapping density function has an uneven distribution. And the mapping is in the case of 30,000 iterations in the [0, 1] interval, and the value density is higher in the [0, 0.1] and [0.9, 1] intervals, reaching 1,900 times. In the interval of [0.1, 0.9], the average value is only 200 times, and the uniform distribution is poor. Similarly, when the Tent map is iterated 30000 times in the [0, 1] interval, the highest value is 439 times, the lowest value is 248 times, and the average is 300 times. The distribution ratio is uniform, which can speed up the search speed. Based on this, the Tent map is selected to generate the turbidity sequence.In the initialization process of Tent chaos,U random integers with a value range of [0, n] are first generated and combined into a chromosome set a1,a2,…,ai,…,aU. Then randomly generate U random numbers with a value range of [0,1], use the Tent mapping expression to chaotically disturb them, generate chaotic variables, and form a chaotic sequence Xij. Then map the chaotic variable to the value range of the optimized variable according to the following formula and round it:(20)gij=1+n−1Xij,i=1,2,…,U;j=1,2,…,k.Obtain the chaotic chromosomeYi=gi1,gi2,…,giki, delete the same genes in them, delete the adjacent 0 when multiple 0s are adjacent, and keep the adjacent one. Calculate the fitness of each chromosome Yi, and select the top TOP M individuals as the initial population. The initial population size in this article is 100.


## 4.2.3. Fitness Function

Fitness function is not only an important basis for genetic operation but also a standard to measure the quality of individuals. The larger the value of the fitness function, the better. However, in practical problems, the objective function is not always nonnegative, sometimes the maximum value and sometimes the minimum value.The greater the fitness of an individual, the better the performance of the individual. However, the objective function of the green and low-carbon cold chain distribution routing optimization problem is the minimum cost. Therefore, the objective function needs to be converted into fitness, and the reciprocal of the objective function can be expressed as the fitness function. And because the optimization problem under the time window constraint is established, the time constraint can be dealt with by adding a penalty function. If the solution corresponding to an individual violates the time constraint, a certain penalty will be given to make it have a smaller fitness, and the probability of being inherited will also become smaller. Therefore, the constructed fitness function is expressed as the inverse of the sum of the objective cost function and the penalty cost of violating the time window, expressed as follows:(21)fi=1zi+ψi,where fi represents the fitness corresponding to individual i, zi represents the objective function value corresponding to individual i, and ψi represents the penalty cost of individual i for violating the time window.


## 4.2.4. Selection Operator

The selection operation determines the probability of being selected according to the size of the fitness value. It is based on the evaluation of the individual’s fitness. The higher the evaluation is, the easier it is to be selected and inherited. The probability of individual selection is calculated as follows:(22)Pk=fk∑i=1nfi.First, calculate the fitness of each body in each generation; sum and then calculate the total value of the individual and the fitness, that is, the probability of the individual’s selection being inherited. Sort all individual selection probabilities. The higher the probability, the better the individual performance, and the higher the ranking, the easier it is to be copied to the next generation.After the selection operation, sort according to the fitness value, and select the firstm−1 chromosomes to perform turbidity perturbation again according to equation (20). Calculate the fitness of each perturbed chromosome, inherit the chromosome individual with the largest fitness value directly to the next generation, and perform the next crossover and mutation operations on the remaining chromosomes. This makes the outstanding individuals after the selection operation chaotic again, increases the diversity of the population, avoids premature maturity, and accelerates the convergence speed of the algorithm.


## 4.2.5. Crossover Operator

The cross-reorganization is carried out using the cyclic crossover rule, and the specific operation process is shown in Figure3.Figure 3
Cyclic crossover operator.Firstly, find a cycle based on the corresponding customer positions in the two individuals of the parent. Then, copy the customer in the cycle position of an individual in the parent to a descendant. Finally, delete the customer in the recurring position of another individual in the parent, and copy the customer in the remaining nonrecurring position directly to the remaining position of the offspring.


## 4.2.6. Mutation Operator

The mutation operation is to change the value of some genes of an individual. The use of mutation operators in GA can improve the algorithm’s ability to solve problems, speed up the generation of optimal solutions, and make the population more diverse and reduce the probability of premature maturity. According to the characteristics of the optimized problem, the use of inversion mutation is also called inverse transformation. Examples are as follows:(1)
Randomly select 2 points from the individual. For example,temp = 589437612; randomly select the 3rd and 6th points; that is, temp=589437612; “∥” represents the crossover sequence;(2)
Perform the reverse transformation operation to obtain the mutated individualtemp1=587349621. The mutation operation is shown in Figure 4.Figure 4
Mutation operation.


## 4.2.7. Termination Criteria

The termination criterion is a condition for judging whether to stop the operation, and the selected termination condition is to reach the preset 2000 generations of evolution. If it reaches the evolutionary algebra, then stop the evolution. In the last-generation population, the distribution path corresponding to the individual with the best fitness value is the optimal solution to the green cold chain logistics distribution routing optimization problem.


## 5. Experiment and Analysis

In order to verify the effectiveness of the proposed algorithm, 25 supermarket stores are selected as delivery customers, the rated load of vehicles is 1T, the fuel type is diesel, the no-load constant velocity fuel consumption is 15.9 L/100 km, and the comprehensive fuel consumption is 24.5 L/100 km, andthe parameter values in the model are shown in Table1.Table 1
Model parameter value.
ParameterValueParameterValuep/(Yuan/t)1200ϑ∗/(L/km)0.382∂1/(Yuan/h)90ϑ0/(L/km)0.168∂2/(Yuan/h)90τ/(kg/L)2.750


### 5.1. Algorithm Running Results

Determine the parameters of the routing optimization GA based on the actual situation: among them, the initial population size is set to 100; the crossover rate and the mutation rate are calculated according to the adaptive probability mentioned in the improved algorithm; the number of iterations is 2000; the number of vehicles is calculated as 3 cars; chromosome length is 15. The running results of the traditional GA algorithm and the proposed improved GA algorithm are shown in Figure5.Figure 5
Running results of traditional GA algorithm and improved GA algorithm.It can be seen from Figure5 that the traditional GA algorithm tends to be stable when iterates about 500 times and reaches the optimal state. Under the same conditions, the improved GA algorithm stabilizes after about 100 iterations and reaches the optimal state. It can be clearly seen that the improved GA algorithm is much faster than the traditional GA algorithm in terms of convergence speed, and the effect of avoiding premature maturity is much better.


### 5.2. Vehicle Route and Distribution Route

Using Matlab7.10.0 (R2010a) programming to solve the green logistics distribution location-routing optimization problem, the resulting distribution path is shown in Figure6.Figure 6
Vehicle distribution route.The three vehicles are the vehicles of the three distribution centers. The specific vehicle running path obtained according to the distribution path is as follows:Vehicle A: 0 ⟶ 7 ⟶ 5 ⟶ 8 ⟶ 11 ⟶ 15 ⟶ 10 ⟶ 19 ⟶ 13 ⟶ 22 ⟶ 25 ⟶ 0Vehicle B: 0 ⟶ 17 ⟶ 14 ⟶ 20 ⟶ 18 ⟶ 23 ⟶ 16 ⟶ 21 ⟶ 2 ⟶ 0Vehicle C: 0 ⟶ 3 ⟶ 12 ⟶ 24 ⟶ 9 ⟶ 4 ⟶ 6 ⟶ 1 ⟶ 0It can be seen that distribution center 0 needs to send out 3 vehicles to serve 25 customers, of which vehicle A serves customers 7, 5, 8, 11, 15, 10, 19, 13, 22, and 25, vehicle B serves customers 17, 14, 20, 18, 23, 16, 21, and 2 services, and vehicle C serves customers 3, 12, 24, 9, 4, 6, and 1; the vehicle returns to the distribution center after completing the distribution service and the total cost of the most distribution process is 12,593 yuan.


### 5.3. Influence Trend of Distribution Penalty Coefficient on Total Cost

In the experiment, it is necessary to analyze the parameters of the proposed algorithm to grasp the influence of the changes of each parameter on the total cost. Due to the perishable nature of fresh food, arriving early during the distribution process will increase the cost of refrigeration due to early refrigeration storage. Delayed arrival will affect the quality of the product and cause penalty costs. Now the distribution penalty coefficient is selected with different values to analyze the importance of the timeliness of fresh food, and the change trend of cost changes is shown in Figure7.Figure 7
Influence of distribution penalty coefficient on the total cost.It can be seen from Figure7 that due to the particularity of fresh food, after the delivery time is exceeded, the food itself will have a penalty cost due to its perishable characteristics. With the increase of the penalty coefficient, the total cost expenditure will also increase. The fluctuation of the penalty coefficient is positively correlated. When the distribution penalty coefficient reaches 4.0, the distribution cost at this time is as high as 13,300 yuan, which is an increase of 707 yuan compared with the normal distribution of 12,593 yuan. Therefore, attention should be paid to the timeliness of the distribution of fresh food, and it is also of practical significance to study the impact of penalty costs on the total cost.


### 5.4. Comparative Analysis of Calculation Results

In the experiment, 25 customers were allocated to 3 distribution centers, a partition distribution mechanism was adopted, and the improved GA algorithm was used to solve the optimal distribution route of each single distribution center. In order to prove the performance of the proposed algorithm, compare it with [14, 17], and the calculation results are shown in Table 2.Table 2
Comparison of calculation results of distribution route.
Distribution centerNumber of usersTraveling distance (km)Total mileage (km)Total cost (yuan)ProposedA1037.92128.9612593B848.19C742.85Reference [17]A732.56134.5513320B1257.08C644.91Reference [14]A631.24137.6613903B538.79C1467.63It can be seen from Table2 that compared with other algorithms, the total mileage and total cost of the proposed algorithm are the lowest, which are 128.96 km and 12,593 yuan, respectively. Because carbon emissions, cargo damage, and other costs are considered, and the routing optimization model is constructed under the multidistribution center mode, the improved GA is used to quickly solve the problem. Therefore, the overall performance is relatively ideal, which can ensure the economy while reducing carbon emissions. Reference [14] realized the global optimal planning of the motion path based on GA but only considered the energy consumption, resulting in vehicle C taking on a large number of distribution tasks, and its overall distribution cost reached 13,903 yuan, an increase of 1,310 yuan compared with the proposed algorithm. Reference [17] proposed a path planning model based on a convolutional neural network, which can effectively improve the efficiency of vehicle path planning. However, the environmental protection and economy of the route have not been considered, so the total mileage is 5.49 km longer than the proposed algorithm.


## 5.1. Algorithm Running Results

Determine the parameters of the routing optimization GA based on the actual situation: among them, the initial population size is set to 100; the crossover rate and the mutation rate are calculated according to the adaptive probability mentioned in the improved algorithm; the number of iterations is 2000; the number of vehicles is calculated as 3 cars; chromosome length is 15. The running results of the traditional GA algorithm and the proposed improved GA algorithm are shown in Figure5.Figure 5
Running results of traditional GA algorithm and improved GA algorithm.It can be seen from Figure5 that the traditional GA algorithm tends to be stable when iterates about 500 times and reaches the optimal state. Under the same conditions, the improved GA algorithm stabilizes after about 100 iterations and reaches the optimal state. It can be clearly seen that the improved GA algorithm is much faster than the traditional GA algorithm in terms of convergence speed, and the effect of avoiding premature maturity is much better.


## 5.2. Vehicle Route and Distribution Route

Using Matlab7.10.0 (R2010a) programming to solve the green logistics distribution location-routing optimization problem, the resulting distribution path is shown in Figure6.Figure 6
Vehicle distribution route.The three vehicles are the vehicles of the three distribution centers. The specific vehicle running path obtained according to the distribution path is as follows:Vehicle A: 0 ⟶ 7 ⟶ 5 ⟶ 8 ⟶ 11 ⟶ 15 ⟶ 10 ⟶ 19 ⟶ 13 ⟶ 22 ⟶ 25 ⟶ 0Vehicle B: 0 ⟶ 17 ⟶ 14 ⟶ 20 ⟶ 18 ⟶ 23 ⟶ 16 ⟶ 21 ⟶ 2 ⟶ 0Vehicle C: 0 ⟶ 3 ⟶ 12 ⟶ 24 ⟶ 9 ⟶ 4 ⟶ 6 ⟶ 1 ⟶ 0It can be seen that distribution center 0 needs to send out 3 vehicles to serve 25 customers, of which vehicle A serves customers 7, 5, 8, 11, 15, 10, 19, 13, 22, and 25, vehicle B serves customers 17, 14, 20, 18, 23, 16, 21, and 2 services, and vehicle C serves customers 3, 12, 24, 9, 4, 6, and 1; the vehicle returns to the distribution center after completing the distribution service and the total cost of the most distribution process is 12,593 yuan.


## 5.3. Influence Trend of Distribution Penalty Coefficient on Total Cost

In the experiment, it is necessary to analyze the parameters of the proposed algorithm to grasp the influence of the changes of each parameter on the total cost. Due to the perishable nature of fresh food, arriving early during the distribution process will increase the cost of refrigeration due to early refrigeration storage. Delayed arrival will affect the quality of the product and cause penalty costs. Now the distribution penalty coefficient is selected with different values to analyze the importance of the timeliness of fresh food, and the change trend of cost changes is shown in Figure7.Figure 7
Influence of distribution penalty coefficient on the total cost.It can be seen from Figure7 that due to the particularity of fresh food, after the delivery time is exceeded, the food itself will have a penalty cost due to its perishable characteristics. With the increase of the penalty coefficient, the total cost expenditure will also increase. The fluctuation of the penalty coefficient is positively correlated. When the distribution penalty coefficient reaches 4.0, the distribution cost at this time is as high as 13,300 yuan, which is an increase of 707 yuan compared with the normal distribution of 12,593 yuan. Therefore, attention should be paid to the timeliness of the distribution of fresh food, and it is also of practical significance to study the impact of penalty costs on the total cost.


## 5.4. Comparative Analysis of Calculation Results

In the experiment, 25 customers were allocated to 3 distribution centers, a partition distribution mechanism was adopted, and the improved GA algorithm was used to solve the optimal distribution route of each single distribution center. In order to prove the performance of the proposed algorithm, compare it with [14, 17], and the calculation results are shown in Table 2.Table 2
Comparison of calculation results of distribution route.
Distribution centerNumber of usersTraveling distance (km)Total mileage (km)Total cost (yuan)ProposedA1037.92128.9612593B848.19C742.85Reference [17]A732.56134.5513320B1257.08C644.91Reference [14]A631.24137.6613903B538.79C1467.63It can be seen from Table2 that compared with other algorithms, the total mileage and total cost of the proposed algorithm are the lowest, which are 128.96 km and 12,593 yuan, respectively. Because carbon emissions, cargo damage, and other costs are considered, and the routing optimization model is constructed under the multidistribution center mode, the improved GA is used to quickly solve the problem. Therefore, the overall performance is relatively ideal, which can ensure the economy while reducing carbon emissions. Reference [14] realized the global optimal planning of the motion path based on GA but only considered the energy consumption, resulting in vehicle C taking on a large number of distribution tasks, and its overall distribution cost reached 13,903 yuan, an increase of 1,310 yuan compared with the proposed algorithm. Reference [17] proposed a path planning model based on a convolutional neural network, which can effectively improve the efficiency of vehicle path planning. However, the environmental protection and economy of the route have not been considered, so the total mileage is 5.49 km longer than the proposed algorithm.


## 6. Conclusion

In the multidistribution center model, optimization goals are designed to minimize costs such as carbon emissions and use improved GA algorithm to solve, finally obtain an economical and environmentally friendly green logistics distribution location optimization path. The proposed algorithm is verified experimentally based on the MATLAB simulation platform. The experimental results show that improved GA has a faster convergence rate than traditional GA, and the distribution penalty coefficient is positively correlated with the total cost. Moreover, the total distance and total cost of the distribution plan obtained by the proposed algorithm are 128.96 km and 12,593 yuan, respectively, which are less than other comparison algorithms. This can provide theoretical support for the distribution routing optimization of logistics enterprises in a low-carbon environment and has a certain reference value.The proposed optimization model sets the corruption rate of fresh agricultural products as a constant and only considers a single type of fresh agricultural products. However, many cold chain products will be distributed together, and the corruption rate is dynamically changing in the actual cold chain logistics distribution. Therefore, how to achieve the distribution optimization goal under the condition of multiple products being distributed at the same time and how the corruption rate is changing need further research.



---

*Source: 7692135-2022-04-12.xml*