#### Abstract

In this study, we investigate a bilevel optimization model for the hazmat transportation problem with lane reservation. The problem lies in selecting lanes to be reserved in the network and planning paths for hazmat transportation tasks. The trade-off among transportation cost, risk, and impact on the normal traffic is considered. By using the traffic flow theory, we quantify the impact on the normal traffic and modify the traditional risk measurement model. The problem is formulated as a multiobjective bilevel programming model involving the selection of reserved lanes for government and planning paths for hazmat carriers. Two hybrid metaheuristic algorithms based on the particle swarm optimization algorithm and the genetic algorithm, respectively, are proposed to solve the bilevel model. Their performance on small-scale instances is compared with exact solutions based on the enumeration method. Finally, the computational results on large-scale instances are compared and sensitivity analysis on the key parameters is presented. The results indicate the following: (1) Both algorithms are effective methods for solving this problem, and the method based on the particle swarm optimization algorithm requires a shorter computation time, whereas the method based on the genetic algorithm shows more advantages in optimality. (2) The bilevel model can effectively reduce the total risk of the hazmat transportation while considering the interests of hazmat carriers and ordinary travellers. (3) The utilization rate of reserved lanes increases with an increasing number of tasks. Nevertheless, once the proportion of hazmat vehicles becomes excessive, the advantage of reducing the risk of the reserved lanes gradually decreases.

#### 1. Introduction

Owing to the rapid development of the chemical industry, the transportation demand for hazmat keeps increasing [1]. Road transportation occupies an important position in the transportation of hazmat and shows a rapid development trend. However, in recent years, many road traffic accidents involving hazmat have brought catastrophic losses to human beings and the environment [2]. There is an average of 36 serious social-impact accidents involving road transportation of hazmat in China annually [3].

To reduce the accident risk, the lane reservation strategy has been employed in hazmat transportation [4–7]. The lane reservation can effectively reduce the risk of hazmat transportation, but it will inevitably affect the driving efficiency and time cost of ordinary vehicles [4]. Thus, it is important to measure the risk and the impact on the normal traffic accurately. Further, the probability of traffic accidents varies with traffic flow [8] and traffic state (e.g., the proportion of hazmat vehicles and the use of lanes) [9]. Also, the impact due to the lane reservation can be improved based on the traffic flow [10]. Therefore, motivated by this dynamic characteristic, this work is an expansion of previous Lane Reservation for Hazmat considering the impact of traffic flow.

Existing researches formulated this problem into a multiobjective model, which ignored the hierarchical and master-slave relationship of decision-making between the government and hazmat carriers [11]. Further studies from the perspective of bilevel optimization are still necessary. Indeed, the lane reservation problem can be regarded as a special form of road network design based on the existing infrastructure. This paper is the modelling integration of hazmat transportation problem and network design problem from policy integration perspective. That can be interpreted as an expansion of the policy integration proposed by [12] considering the exclusive right-of-way of hazmat vehicles. Additionally, there is still dearth of efficient heuristic algorithms for large-scale problems [6].

Based on these gaps, this study focuses on a bilevel optimization model for the hazmat transportation problem with consideration of the lane reservation strategy. The contributions are as follows: (1) We propose a bilevel programming model for hazmat transportation problem via lane reservation. (2) Based on the traffic flow density, vehicle speed, and other specific traffic information, we improve the quantification of the transportation risk and the impact on the normal traffic, which makes the model more practical. (3) Two hybrid metaheuristic algorithms are proposed and compared with the accurate algorithm based on enumeration method to verify the validity of the model and algorithm. (4) The effects of the number of hazmat transportation tasks and the proportion of hazmat vehicles are analyzed to provide some management insights.

The rest of the paper is organized as follows. A literature review from two aspects, the hazmat transportation and lane reservation problem, is presented in Section 2. The problem description and formulation are described in Section 3. Two heuristic algorithms are proposed in Section 4. In Section 5, computational experiments are conducted to demonstrate the effectiveness of model and algorithms. Finally, Section 6 presents conclusions and potential directions for future research.

#### 2. Literature Review

##### 2.1. Hazmat Transportation

The hazmat transportation problem is concerned with finding routes with minimum cost and minimum risk [13]. Our problem is considered as a transportation network design problem. Its mathematical programming formulation generally relates to regulatory policies on shipments with hazmat, aiming at minimizing the total network risk [12].

The most important way is to restrict hazmat vehicles to travel on certain sections. Kara and Verter [14] first proposed the problem of designing a network by imposing curfews. Then He et al. [11] investigated the hazmat transportation problem with traffic restriction addressing also the concept of risk equity and provided a heuristic algorithm to always find a stable solution. Taslimi et al. [15] established a bilevel network design model considering the response time of accidents and risk equality. Esfandeh et al. [16] extended the network design problem to account for time-dependent road closure.

Additionally, other traffic management policies have been integrated with network design problem. Marcotte et al. [17] proposed the toll policies to guide hazmat vehicles to route according to their economic and time preferences, thereby deterring hazmat vehicles from using certain road segments. Esfandeh et al. [18] studied the dual toll solution for hazmat vehicles and ordinary vehicles to regulate hazmat transportation. Wang [19] introduced the speed limit strategy for hazmat transportation. López-Ramos et al. [12] presented an integrated model and specialized local search for the road network pricing and the design of hazmat vehicles.

##### 2.2. Lane Reservation Problem

The lane reservation is to select the lanes and/or the time intervals to reserve for only special tasks to satisfy the condition of time or safety. It was initially applied to the strategic planning of transportation tasks for large-scale events [20]. Wu et al. [21] first studied the lane reservation problem in time-constrained transportation networks as a mathematical model and obtained near-optimal solutions. Fang et al. [22] investigated a capacitated lane reservation problem by considering the residual capacity of a lane and consequently developed an exact algorithm based on cut-and-solve approach. Fang et al. [23] extended an optimal lane reservation problem with a dynamic link travel time which is closer to a real traffic situation. The most outstanding application of lane reservation is the exclusive bus lanes design [24]. In recent years, it has been gradually applied to automated truck freight transportation [25].

Several studies have been carried out on lane reservation for hazmat. Zhou et al. [4] first introduced the lane reservation strategy into hazmat transportation. A multiobjective model is developed to minimize the total transportation risk and the impact on the normal traffic. Then they proposed a *ε*-constraint and fuzzy logic-based optimization for this problem [5]. Chu et al. [6] investigated the time-dependent lane reservation for hazmat and solved it via a cut-and-solve-based *ε*-constraint method. Then a new cut-and-solve and cutting plane combined method is proposed to reduce the computational time [7].

#### 3. Mathematical Model

##### 3.1. Problem Description

Let be a hazmat transportation network, where is the set of nodes and is the set of arcs; arc denotes a road segment from node to node . Hazmat transportation tasks are described by sets of origins () and destinations ().

The government has a leading position concerning the hazmat transportation risk, whereas hazmat carriers attempt to minimize their cost. We design a bilevel hazmat lane reservation problem (BHLRP) as depicted in Figure 1. The upper level (UP) means the government selects lanes to be reserved considering the interests of hazmat carriers and ordinary travellers. That is done to minimize the transportation risk and the impact on the normal traffic. The lower level (LP) is the optimal path with the lowest cost for each task.

##### 3.2. Assumption and Notation

To better define the problem, we present some assumptions and notations here.

###### 3.2.1. Assumptions

(1)Each road section has at least two lanes so that at least one general lane is in use after lane reservation(2)When a lane is reserved in the section, hazmat vehicles passing through a road must use the reserved lane(3)Vehicles are traveling in each section at a constant speed(4)Potential accidents involving hazmat independently occur, and the probability of hazmat accidents is affected by traffic flow density(5)There are only two types of vehicles in the road network, namely, hazmat vehicles and ordinary vehicles

Decision variables:

##### 3.3. Impact on Normal Traffic

With the development of intelligent transportation network, traffic situations such as speed and traffic density can be obtained by analyzing huge amounts of real traffic data. This provides a new way to further quantify the impact on the normal traffic. We improve this indicator based on the traffic flow by referring to [26].

Before implementing the lane reservation strategy, vehicles are evenly distributed on each lane of arc and travel at the same speed of . After lane reservation, the traffic d with travel speed seen in equation (2) is derived from the speed-density relationship of traffic flow theory. is the optimal speed with the maximum flow. Thus, the impact on the normal traffic in arc , , can be defined by equation (3).

##### 3.4. Risk Assessment

The risk assessment model used in our previous work, which considers risk as the product of the probability of hazmat accident and its consequence, ignores the impact of traffic flow density on the probability of accidents.

Taking traffic flow density into account, we improve the method to measure the risk of hazmat vehicles passing reserved lanes. Before lane reservation, hazmat vehicles travel on each lane of arc with traffic density of and accident probability of . After lane reservation, hazmat vehicles travel on the reserved lane and the probability of hazmat accident can be defined as

##### 3.5. Mathematical Model

In this subsection, we formulate an integer linear bilevel programming model for BHLRP; that is, the objective functions and the constraints are linear. All three decision variables are binary. Moreover, this formulation is under the partial cooperation assumption [27], also called an optimistic BHLRP. This implies that the government assumes that hazmat vehicles’ choice is always the most favorable one when there are multiple costly equivalent alternative travel routes. The model is presented as follows.

The UP of BHLRP:

The LP of BHLRP:

Equations (5)–(8) indicate the UP of the model and equations (9)–(13) show the LP of the model. The objective function (5) is to minimize the transportation risk. The objective function (6) is formulated to minimize the impact on normal traffic caused by the lane reservation strategy. Constraint (7) ensures the total risk should be reduced by *α* after lane reservation. Constraint (8) shows the 0-1 restriction on the upper decision variable.

The objective function (9) seeks to minimize the total time cost of hazmat transportation. Constraint (10) guarantees that hazmat tasks should be transported from origination nodes to destination nodes and the flow balance at each node. Constraint (11) means that a reserved lane may be used for task only when the task’s path passes the arc . Constraint (12) ensures that a reserved lane may be used by task only when the lane is reserved. Constraint (13) represents 0-1 restriction on the lower decision variables.

#### 4. Solution Algorithm

This paper formulates a bilevel optimization model for BHLRP. In addition to the complexity of the routing problem, it is understood that the bilevel model is NP-hard [28]. Several exact methods depend on some characteristics of the solution space, such as the vertex enumeration [29] and the KKT method [11, 14]. These can only solve partial bilevel programming problems [30]. Furthermore, the heuristic algorithm, such as the particle swarm optimization (PSO) [31, 32] and the genetic algorithm (GA) [33, 34], is applied to solving bilevel problems because of its better global search capability and its lower dependence on the form of the objective function.

In the bilevel model, we will consider the lane reservation and routing problems simultaneously. The UP has decision priority and can predict the response of the LP to each decision. The LP is a multisource shortest path problem that obeys the shortest path allocation principle. Therefore, we propose a three-stage solution to resolve this problem, presented in Table 1. According to the route of tasks, a 0-1 matrix () is introduced in step 1 to represent candidate sections to be reserved. Because there is a task passing through arc , that is, , it is a candidate road section with a reserved lane (=1).

To obtain an approximate optimal solution in a reasonable time, based on the three-stage solution, we design two hybrid heuristic algorithms named EX-PSO-EX and EX-GA-EX, combining PSO and GA with the shortest path algorithm, respectively, to solve BHLRP. To verify the effectiveness of the proposed algorithms, there is a benchmark algorithm combining explicit complete enumeration with the shortest path algorithm, named EX-EX-EX. This method is only suitable for solving small-scale examples. The optimal solution obtained by it can be compared to the approximate optimal solution obtained by heuristic algorithms to evaluate their performance. Table 2 illustrates the general structure of the three approaches.

##### 4.1. Solution Code

In BHLRP, the decision variable of UP is binary; thus, a type of layout of reserved lanes can be represented as a binary matrix of two dimension . The position equals 1 if there is a reserved lane on arc and 0 otherwise. Figure 2(a) describes a solution to the lane reservation in a 10-node network. As seen in this reservation plan, arcs (4, 5), (4, 6), (6, 4), (7, 8), (8, 9), and (8, 10) are reserved.

**(a)**

**(b)**

The decision variable of LP is notated as a 0-1 matrix of three dimensions . Each page denotes the path of task . The position equals 1 if task passes arc but 0 otherwise. The path of a task from node 10 to node 4 is displayed in Figure 2(b). The path is 10-8-6-4. The decision variable of LP , obtained by the dot product of and , is also expressed in a 0-1 matrix of three dimensions .

##### 4.2. Fitness Calculation

First, objective functions (6), (6), and the left part of constraint (7) are marked as , and , respectively. As the dimensions of and are not uniform, after maximum-minimum normalization, we record their linear weighted sum as the new objective function. Then is dealt with via the penalty function. Hence, the fitness function can be expressed aswhere is a sufficiently large positive number to indicate the penalty factor.

##### 4.3. EX-PSO-EX

In the EX-PSO-EX algorithm, a Sigmoid function is introduced to deal with the velocity of particles when determining the reserved lanes via PSO. Thus, we can ensure the position value of each particle is 0 or 1. Further, the position updated formula is improved as follows:

The Sigmoid function implies that the probability of value in the next step is 1 by mapping particle velocity to interval. Moreover, ensures that the particles update only at the position of candidate sections to be reserved. It ensures the updated particles are still relatively optimal feasible solutions.

Additionally, adaptive weight and asynchronous learning factors are used to balance the global search ability and local improvement ability of the PSO algorithm. The updated formulas of the weight and the learning factor are shown in (16) and (17), respectively.

Let and , respectively, represent the maximum and minimum values of . is the current fitness value of the particle. Moreover, and indicate the average fitness value and the minimum fitness value of all current particles, respectively.where and are the initial and final values of the iteration of , respectively, is the current number of iterations, and denotes the maximum number of iterations. We refer to literature [35] to obtain their values: , , , and . Figure 3 shows the pseudocode to solve BHLRP by using the EX-PSO-EX method.

##### 4.4. EX-GA-EX

By using GA to determine lanes to be reserved in step 2, we propose the second hybrid metaheuristic method, EX-GA-EX. The pseudocode is shown in Figure 4. In EX-PSO-EX, binary coding is adopted to code each chromosome into a matrix of , representing a solution for lane reservation.

To search for the solution space effectively, we consider the genetic operators, namely, the ranking selection operator, uniform cross-block operator, and uniform mutation operator, as illustrated in Figure 5. To ensure the offspring after crossing and mutation are still a relatively better feasible solution, both the intersection block and mutation location are generated in the candidate segments of lane reservation.

**(a)**

**(b)**

##### 4.5. Benchmark Algorithm: EX-EX-EX

In the EX-EX-EX method, we enumerate the solutions of the lane reservation in step 2 completely based on the candidate segments obtained in step 1. To some extent, it minimizes the enumeration space and enumeration quantity. The pseudocode of EX-EX-EX is given in Figure 6.

#### 5. Computational Results

In this section, computational results are presented to verify the effectiveness of the above model and algorithms and to assess the impact of the key factors. All the proposed hybrid algorithms are encoded with MATLAB 2018 software and all experiments are performed on a server with 3.00 GHz speed and 8 GB RAM.

Based on random network topology generator [36], the road network is randomly generated. This approach is presented in Table 3. Then 60 instances of different scales are created ranging from 8 to 100 nodes and from 5 to 30 tasks and different proportions of hazmat vehicles of 0.1-0.2. To ensure the best performance of PSO and GA, combining with the Taguchi Method [29], the parameters for different scale examples are calculated through a large amount of presolving, as shown in Tables 4 and 5.

##### 5.1. Comparison of Algorithm Performance

After generating 60 random instances, instances in small scale () were solved based on the EX-EX-EX approach. Then the EX-PSO-EX and EX-GA-EX approaches were applied to solve all instances. Each instance is solved with different algorithms 15 times and their averages are recorded.

The results for small-scale instances are listed in Table 6. Columns 2–4 contain the sizes of the instances. The next column shows the number of candidate sections to be reserved (). Columns 6-7 show the upper objective function value () of the optimal solution () and CPU time obtained by EX-EX-EX. Columns 8–19 show the average value of (Average Solution), the percentage deviation of the solution’s mean with respect to (Avg), the minimum of (Best Solution), the percentage deviation of the solution’s minimum with respect to (Best), the standard deviation of (Standard Deviation), and the CPU time obtained by EX-PSO-EX and EX-GA-EX, respectively.

Note that the CPU time required by EX-EX-EX is not only affected by , , and but also directly related to , exhibiting exponential growth in small-scale examples. Set 15 with only 15 nodes, 78 sections, and 17 candidate sections has been running for up to 1381.517 s. However, with EX-PSO-EX and EX-GA-EX, all sets in Table 6 can obtain their effective approximate optimal solution with CPU time of less than 2 s, which is much lower than EX-EX-EX. Although EX-GA-EX takes a while longer than EX-PSO-EX in most of the examples, the optimal solution of 13 examples is found in each of the 15 runs in Table 6. From the other 6 indicators, such as Average Solution and Avg, the performance of EX-GA-EX is better than that of EX-PSO-EX.

To further test the performance of EX-PSO-EX and EX-GA-EX, we use them to solve 15 larger examples. The CPU time to solve them with EX-EX-EX is more than 3600 s. Their Average Solution, Best Solution, Standard Deviation, and the CPU time of each instance are shown in Table 7. Additionally, Figure 7 shows the time for solving the 3 methods proposed. Furthermore, the Standard Deviations of sets 1–30 obtained by EX-PSO-EX and EX-GA-EX are illustrated in Figure 8.

In terms of average CPU time, it takes 85.214 s to solve set 30 with EX-GA-EX, while EX-GA-EX consumes about half the time taken by EX-PSO-EX. However, the performance of EX-GA-EX in the optimal and average solution quality is better. The maximum standard deviation is only 0.014. The results show that, even in larger instances, the performance of EX-PSO-EX and EX-GA-EX is robust, which can be used as an effective method to solve the problem of lane selection for BHLRP.

##### 5.2. Sensitivity Analysis

To make a trade-off between the effect of risk reduction and the impact on the normal traffic, it is necessary to further evaluate the effectiveness of lane reservation solutions with different numbers of tasks and hazmat vehicle ratios. Two scenarios are set up and each scenario tests several instance sets.

Considering the good performance of EX-GA-EX in different scales, we use it to solve each set 15 times and record and use the average value obtained. The evaluation results are shown in Tables 8–11 and Figures 9-10.

###### 5.2.1. Effect of the Number of Tasks ()

The number of tasks reflects the size of the examples. Table 8 shows average computational results of the scenario where varies from 5 to 30 in three networks, 20 × 126, 30 × 286, and 40 × 526. The results include the number of reserved lanes, the risk, the decreased risk compared to the risk before lane reservation (), the relative change in risk, the impact on normal traffic (), and . Firstly, after lane reservation, the value of is much smaller than 1, the value of before lane reservation, and the relative reduction of the risk exceeds 50%. This indicates that the effect of lane reservation in reducing the risk of hazmat transportation is very significant. Secondly, with the increase in , reserved lanes need to be set up on more sections. At this time, the risk, , and show an obvious increasing trend; however, the effect of reserved lanes to reduce the risk is increasingly obvious. The increase in risk is mainly due to the increase in the total risk with the increase in . It is not difficult to imagine that if there are no reserved lanes, then total risk will increase more; the increase in and is due to the increase in the number of road sections with reserved roads.

To evaluate the effectiveness of reserved lanes with varied |W|, Table 9 shows the ratio of and () under three different networks and their average ratio. Figure 9 shows the trend of the average ratio. The results show that, with the increase in , there is a trend of increasing after decreasing. When is small, the decrease shows that reducing the same risk will have a greater impact on the normal traffic. This is mainly because fewer sections can be shared between tasks so that additional reserved lanes need to be set up to meet the risk reduction needs. Then the number of affected sections will increase. When increases to 15 or more, rises instead. It is not possible to set up new reserved lanes. There will be more of same road sections between tasks; that is, the utilization rate will be improved. Therefore, the utilization rate of reserved lanes is very important to ensure the effectiveness of the lane reservation.

###### 5.2.2. Effect of the Proportion of Hazmat Vehicles ()

The proportion of hazmat vehicles will affect the traffic flow density on these two kinds of lanes after lane reservation. It will have an impact on the risk and the normal traffic. To evaluate the effectiveness of reserved lanes with varied , Table 10 shows the average computational results of the scenario where varies from 0.1 to 0.2 in three networks, 20 × 126, 30 × 286, and 40 × 526. With the increase in , and the risk increase, while and the reduced risk per lane () decrease, which reflects the decline in the effect of reserved lanes on risk reduction. The main reason is that more hazmat vehicles can travel on reserved lanes. The traffic flow density on reserved lanes will increase relatively. It will increase the probability of accident to some extent. Reduction of can be attributed to the fact that fewer ordinary vehicles travel at a slower speed. Additionally, the traffic flow density on GP lanes decreases, which will reduce the impact on the traffic time on GP lanes.

The result shows that the number of reserved lanes does not necessarily increase with the proportion of hazmat vehicles. This indicates that more sections are affected. Thus, the time for ordinary vehicles to travel through these lanes will be longer. Additionally, if is increased to a certain extent, the advantage of reserved lanes to reduce risk will be no longer obvious. This also reflects the validity of the model to a certain extent.

Table 11 shows under three different networks with varied and the average ratio. Figure 10 shows the trend of the average ratio. With the increase in *β*, decreases from 4.525 to 3.523. Although the utilization of reserved lanes will improve, the optimization effect becomes decreasingly significant. When the proportion is sufficiently large to have a negative impact on the free flow on reserved lanes, the advantage of lane reservation to risk reduction will gradually decrease. Especially in some road sections where , the traffic flow density on the reserved lane after lane reservation will be greater than that on each lane before lane reservation. This may signify that it may not be sufficient to set up only one reserved lane in these sections.

#### 6. Conclusion

In this study, the problem of hazmat transportation with lane reservation is studied. A bilevel multiobjective model is established for the problem. The measure of the risk and the impact on the normal traffic is improved by incorporating specific traffic information, such as traffic flow. Subsequently we design two hybrid metaheuristic algorithms and an exact algorithm to solve the model. To validate the model and algorithms, 60 examples are created and the sensitivity analysis of two key factors (the number of transportation tasks and the proportion of hazmat vehicles) is conducted. The main conclusions are as follows:(1)The bilevel optimization model for the hazmat transportation problem with the lane reservation can significantly reduce the risk of the hazmat transportation while considering the interests of hazmat enterprises and ordinary travellers.(2)The proposed hybrid metaheuristic algorithms can obtain effective approximate optimal solutions while greatly reducing the solution time. EX-PSO-EX is superior to EX-GA-EX in CPU time. EX-GA-EX has more advantages in the average and optimal solution quality of the approximate optimal solution.(3)The utilization of reserved lanes will affect the effectiveness of risk reduction. To a certain extent, the increase in transportation tasks will improve the utilization of reserved lanes. When the proportion of hazmat vehicles is enough to have a negative impact on the free flow on the reserved lane, the advantage of reducing the risk of the reserved lane will gradually decrease.(4)In sections where hazmat vehicles are particularly dense, , it may not be enough to set up only one reserved lane. The government should consider the road network structure, the vehicle composition ratio, and other relevant traffic data when selecting reserved lanes to further weigh the effect of risk reduction and the impact on normal traffic.

The model we propose does not consider the impact of time-varying conditions and other uncertainties on transportation risk, neither does it consider the impact of social vehicles on hazmat vehicles. Therefore, making the model closer to the actual situation is a further research direction. Another future work may be to consider the integration of lane reservation and other regulation policies. Additionally, in order to effectively reduce the risk and improve the utilization of existing transportation infrastructure of urban areas, a dynamic lane reservation for hazmat transportation may be a more applicable choice based on the actual traffic state and the regulation of traffic state information.

#### Notations

*Sets*

: | Set of hazmat tasks, |

: | Set of origination nodes of tasks, |

: | Set of destination nodes of tasks, . |

*Parameters*

: | Total number of lanes on arc |

: | Distance of arc |

: | Traffic density in the arc |

: | Travel speed in arc without reserved lanes |

: | Travel speed on reserved lanes of arc after lane reservation |

: | Travel speed on general lanes of arc after lane reservation |

: | Travel time in arc without reserved lanes |

: | Travel time on reserved lanes of arc after lane reservation |

: | Travel time on general lanes of arc after lane reservation |

: | Impact on normal traffic due to lane reservation on arc |

: | Population exposure if an accident happens on arc |

: | Accident probability of hazmat happening on the GP lanes of arc |

: | Accident probability of hazmat happening on the reserved lane of arc |

: | The threshold of risk reduction after lane reservation |

: | The proportion of hazmat vehicles in total motor vehicles. |

#### Data Availability

The data used to support the findings of this study are available from the corresponding author upon request.

#### Conflicts of Interest

The authors declare that they have no conflicts of interest.

#### Acknowledgments

This research was supported by the National Natural Science Foundation of China (NSFC) (no. 71901037) and the Fundamental Research Funds for the Central Universities (no. 300102239102).