Collaboration and Resource Sharing in the Multidepot Multiperiod Vehicle Routing Problem with Pickups and Deliveries

: In this work, a multidepot multiperiod vehicle routing problem with pickups and deliveries (MDPVRPPD) is solved by optimizing logistics networks with collaboration and resource sharing among logistics service providers. The optimal solution can satisfy customer demands with periodic time characteristics and incorporate pickup and delivery services with maximum resource utilization. A collaborative mechanism is developed to rearrange both the open and closed vehicle routes among multiple pickup and delivery centers with improved transportation efficiency and reduced operational costs. The effects of resource sharing strategies combining customer information sharing, facility service sharing, and vehicle sharing are investigated across multiple service periods to maximize resource utilization and refine the resource configuration. A multiobjective optimization model is developed to formulate the MDPVRPPD so that the minimum total operational costs, waiting time, and the number of vehicles are obtained. A hybrid heuristic algorithm incorporating a 3D clustering and an improved multiobjective particle swarm optimization (IMOPSO) algorithm is introduced to solve the MDPVRPPD and find Pareto optimal solutions. The proposed hybrid heuristic algorithm is based on a selective exchange mechanism that enhances local and global searching capabilities. Results demonstrate that the proposed IMOPSO outperforms other existing algorithms. We also study profit allocation issues to quantify the stability and sustainability of long-term collaboration among logistics participants, using the minimum costs remaining savings method. The proposed model and solution methods are validated by conducting an empirical study of a real system in Chongqing City, China. This study contributes to the development of efficient urban logistics distribution systems, and facilitates the expansion of intelligent and sustainable supply chains. and J.F.; resources, Y.W.; data curation, Y.W. and Q.L.; writing—original draft preparation, Y.W. and Q.L.; writing—review and editing, Y.W., Q.L., and X.G.; visualization, Q.L. and Y.L..; supervision, Y.W., J.F., and H.W..; project administration, Y.W., J.F., and H.W.; funding acquisition, Y.W. All authors have read and agreed to the published version of the manuscript.


Introduction
Reverse logistics has been widely adopted in the sustainable development of modern logistics, where the pickup and delivery systems are both incorporated to reduce the transportation cost and meet the regulations on environmental protection [1,2]. In previous studies, the vehicle routing problem with pickup and delivery (VRPPD) has received widespread attention on account of its roles in the practical application of reverse logistics [3,4]. Grocery store chains, for example, require services combining delivery and pickup in a single trip. The delivery of bottled drinks and pickup of empty bottles are operated by the same vehicle on a single trip. Many companies have begun controlling the complete life cycle of their products; in particular, organizers have taken it upon themselves to mitigate the environmental impact a product has during its lifetime from being produced to being consumed/disposed of/recycled. In those cases, vehicles are more effectively utilized when performing commodity distribution to customers (i.e., forward logistics) [5] and returning products to the depot (i.e., reverse logistics) [6]. Therefore, how to coordinate and optimize pickup and delivery services in a sustainable logistics system is an issue of practical significance.
However, a challenge in addressing this issue is that the delivery and pickup schedules could generate additional operational costs if they are fixed while the product demand from customers changes and differs on any given service day [7]. The periodic vehicle routing problem (PVRP) and multiperiod vehicle routing problem (MPVRP) are designed to find optimal routes in multiple time periods (i.e., multiple days) with minimum travel distances [8][9][10]. Given the extension of the scale of logistics systems, researchers, and companies have studied periodic routes from one depot or among multiple depots [11,12], and open-closed mixed vehicle routing problems have been devised to formulate trips resembling real-world applications in multidepot logistics networks [13,14]. The extended logistics network structures composed of multiple pickup and delivery depots, however, are usually not optimal in terms of operational efficiency and resource configuration in previous studies. A multidepot multiperiod vehicle routing problem with pickup and delivery (MDPVRPPD) is proposed in this study, not only to design the periodic vehicle routes over a planning horizon for customers, but also to incorporate pickup and delivery services with open and closed routes for network optimization.
From the perspective of sustainable transportation and logistics, collaboration is key to improve travel efficiency, reduce operational costs, and achieve sustainable development [15,16]. The concept of collaboration is widely used in city logistics, including various industries and manufacturers [17,18], urban good transport activities [19], and e-commerce logistics trading systems [20]. Multiple sharing mechanisms including information [21], facilities [22], vehicles [23] and resource [24] sharing have been widely used as an effective method to optimize such logistics operations. The cost allocation mechanism has been studied over the last several decades to allocate or share the profits or cost savings of each participant and maintain collaborative stability [20,25]. Vehicle sharing has been used in the VRPPD as a green public transportation mode to share vehicles for daily travel with high vehicle fill rates [26][27][28]. Consequently, collaboration plays an important role in the transportation system optimization to achieve sustainable logistics development, and resource sharing can be utilized in collaborative logistics pickup and delivery networks to ensure collaborative stability with increased resource utilization.
Collaboration and resource sharing can be mathematically modeled for a multidepot multiperiod logistics network with pickup and delivery (MDPLNPD) optimization. In the present study, a collaborative mechanism and resource sharing strategy are constructed to determine optimal open-closed mixed vehicle routes, coordinate pickup and delivery services with improved operational efficiency and reasonable resource configuration, and, ultimately, solve the MDPVRPPD. A multiobjective integer programming model is developed to formulate the collaborative MDPLNPD and achieve the multiple objectives of minimum total operational costs, number of vehicles, and waiting times (WTs). A hybrid heuristic algorithm combining a 3D clustering method and an improved multiobjective particle swarm optimization (IMOPSO) algorithm is proposed to find optimized vehicle routes and Pareto optimal solutions and develop high-quality optimal solutions for this type of problem. A real-life logistics network with pickups and deliveries is then analyzed to study the collaborative mechanism and resource sharing strategies given multiple facilities and service periods, and demonstrate the validity and applicability of the proposed solution methodology. Optimization of practical networks with minimum operational costs and maximum resource utilization can provide strategic advice for logistics enterprises and promote the sustainable development of urban logistics distribution. utilization and decrease traffic congestion; in this case, vehicles can be used to perform pickup and delivery trips among multiple depots with optimal transportation efficiency [27].
Multiple heuristic or metaheuristic algorithms have been developed by researchers in the VRP to obtain high-quality or optimal solutions with minimal computation times [44]. Goksal et al. [45] introduced a hybrid algorithm by combining PSO and variable neighborhood descent, a well-known metaheuristic solution approach, to address the VRP with simultaneous pickup and delivery (VRPSPD). Yanik et al. [29] developed a hybrid metaheuristic approach combining a genetic algorithm (GA) and a modified savings algorithm to address the VRPPD with time windows. Kalayci and Kaya [11] introduced a hybrid metaheuristic approach composed of ant colony system (ACS) and variable neighborhood search to optimize pickup and delivery services in an optimized network with optimal vehicle schedules. Li et al. [46] developed a mixed-integer linear programming model to formulate the pickup and delivery problem and presented an adaptive large neighborhood search to solve the optimization problem. An integer linear programming model-based local search procedure was proposed by Montero et al. [3] to address the VRPPD. Guo et al. [30] established a mixed-integer programming model to formulate the VRPPD with multiple objectives and obtained higher efficiency and lower complexity compared with those achieved by GA.
A series of vehicle routes should be found to satisfy all the customer demands and service time windows in the VRP with multiple periods and a minimum number of vehicles and travel distances [7,47]. An improved ɛ-constraint method was studied to balance trade-offs in the multiobjective MDPVRP and optimize the interdepot routes of a sustainable reverse logistics system [32]. A hybrid genetic search with adaptive diversity control was proposed by Vidal et al. [47] to solve the MDPVRP considering its computational efficiency and solution quality. A three-stage algorithm was developed by Luo et al. [7] to find high-quality initial solutions for the MPVRP through a decomposition algorithm, a tree search algorithm, and a tabu search (TS) algorithm. Azad et al. [36] introduced a heuristic initialized stochastic memetic algorithm to rebalance exploration and exploitation, and avoid premature convergence in the MDPVRP optimization. Hernandez et al. [33] used a unified TS that evaluates each tactical routing plan in the multiperiod logistics network and obtains the minimum cost under the constraints of customer demands and service times. Estrada-Moreno et al. [34] introduced a two-stage approach consisting of iterated local search and biased-randomization technologies to solve the MPVRP with optimal travel routes.
Despite the wealth of knowledge that has since been obtained, however, the previous research related to the MDPVRPPD presents the following limitations. (1) The collaborative mechanism is rarely considered in logistics network optimization with pickups and deliveries. (2) Resource sharing strategies have been inadequately discussed in multidepot and multiperiod networks. (3) A valid model for optimizing both open-closed mixed vehicle routes in collaborative logistics networks with multiple depots and service periods is currently lacking. (4) State-of-the-art evolutionary solution methodologies have limited applicability in solving the MDPLNPD through collaboration and resource sharing.
The contributions of this study are as follows. (1) A collaborative mechanism involving facilities providing pickup and delivery services is developed to improve the transportation efficiency in the MDPLNPD. (2) Resource sharing strategies are introduced to the MDPVRPPD to maximize resource utilization and facilitate sustainable logistics development. (3) A multiobjective integer programming mathematical model is developed to formulate the open-closed mixed routing problem in the collaborative network with pickups and deliveries. (4) A case study is conducted to demonstrate the effectiveness of the proposed hybrid heuristic algorithm in solving the MDPVRPPD with collaboration and resource sharing.

Problem Statement
In contrast to traditional logistics networks with pickups and deliveries, the MDPLNPD incorporates open and closed routes and aims to optimize the network through effective collaboration and resource sharing among multiple facilities and service periods. Multiple service periods can be determined on the basis of the different time windows of facilities and customers.
Practical network optimization considers transportation among multiple logistics facilities in the first echelon, routes from delivery centers (DCs) or pickup centers (PCs) to customers in the second echelon, and the logistics resource utilization among multiple service periods. The comparison of the noncollaborative and collaborative MDPLNPD is shown in Figure 1. Here, rather than their usual initial independent operation, trucks perform transportation among collaborative facilities, and the routes plied by vehicles for pickups and deliveries are optimized with minimal crisscross and longdistance travels, and the number on each route expresses the travel time by trucks or vehicles. A comparison of networks before and after optimization is shown, and the number indicated in each route represents the travel time of trucks or vehicles. In Figure 1a, the delivery and pickup services are independent among the DCs and PCs in the noncollaborative network. The vehicles used in the pickup and delivery routes are independent, resulting in increased crisscross transportation and reduction in vehicle resource utilization. Unreasonable vehicle scheduling means customers are served beyond their time windows. For instance, C2, C14, C17, and C23 are visited beyond their time windows. Therefore, the noncollaborative logistics network should be optimized with improved transportation efficiency and resource utilization. In Figure 1b, collaboration is established among the DCs and PCs, and customers are reassigned to their adjacent facilities and service periods. Pickup and delivery services are incorporated with open-closed mixed routes. Unreasonable transportation is decreased, and customers are visited within their time windows. Vehicles can be shared among collaborative DCs and PCs in the first and second periods. For example, the closed route (e.g., DC2→C7→C1→C2→ DC2) from DC2 in the first period can be shared with the closed route (e.g., PC1→C22→C21→C20→ PC1) from PC1 in the second period. The open delivery route (e.g., DC1→C8→C14) departure from DC1 in the first period can be shared with the open pickup route (e.g., C17→C16→C15→PC1) back to PC1 in the second period. Therefore, collaboration and resource sharing can be adopted in the network optimization to achieve transportation efficiency and a reasonable resource configuration. DC2 PC1 DC1 [5,7]   (a) Non-collabor ative network before optimization DC2 PC1 DC1 [5,7] [ Time windows of customers i n the 1st and 2nd periods C1 [8,10] Assume that the per time unit cost of pickup or delivery by vehicles is $20, the per unit vehicle maintenance cost is $10, the per time unit penalty cost is $5, and the per unit facility government subsidy (GS) is $60. A comparison of the pickup time (PT), delivery time (DT), WT, number of vehicles, pickup cost (PC), delivery cost (DC), maintenance cost (MC), penalty cost, GS, and total operational cost (TC) in noncollaborative and collaborative MDPLNPD is shown in Table 1.  Table 1 shows that the optimization of the MDPLNPD can reduce the total operational cost and improve the operational efficiency of the logistics system. The total logistics operational cost decreases from $1925 to $1470 in the two service periods considered, the pickup and delivery costs reduce by 31.7% and 34.0%, and the maintenance cost declines from $130 to $50 on the condition of the collaboration and resource sharing among multiple facilities. No waiting time is found in the collaborative network, and the penalty cost is minimized after network optimization. The number of vehicles is also reduced from 13 to 5 on account of vehicle sharing among multiple facilities and service periods. Therefore, effective collaboration and resource sharing can improve operational efficiency and maximize resource utilization.

Model Formulation
The mathematical model of collaboration and resource sharing for MDPVRPPD optimization is described in this section. The related definitions and parameters used to formulate the MDPVRPPD are presented in Tables A1 and A2 in Appendix A.
A multiobjective mathematical programming model is constructed to achieve MDPLNPD optimization through effective collaboration and resource sharing. Three objectives, namely, minimum logistics operational cost T, waiting time Z, and number of vehicles H, are established. In Equation (1), the total logistics operational cost TC includes the operational cost TC1 among the logistics facilities in the first echelon and the operational cost TC2 among the logistics facilities and customers in the second echelon. Equation (2) is proposed to minimize the total waiting times Z1 for delivery and Z2 for pickup in the MDPLNPD. Equation (3) determines the minimum number of vehicles on the basis of the maximum number of shared vehicles within multiple service periods. 1 2 is the delivery and pickup cost generated by vehicles traveling between facilities and customers.
is the penalty cost of earliness or delay when vehicles arrive at the delivery customers.
presents the penalty cost of earliness or delay when vehicles arrive at the pickup customers.
as the waiting time of vehicles at the delivery customers.
Z2 shown in Equation (7) presents  as the waiting time of vehicles at the pickup customers.
Constraints in the first echelon: , , , , , , , , Constraint (8) denotes that each truck can only depart from one facility. Constraints (9) and (10) ensure that each DC can be visited by a single truck from another DC and each PC can be visited by a single truck from another PC within the kth period. Constraints (11) and (12) imply that each truck leaves to serve another DC after serving one DC and each truck leaves to serve another PC after serving one PC. Constraints (13) and (14) can be used to eliminate subtours on every route departure from DCs or PCs. Constraints (15) and (16) guarantee that each truck has sufficient transport capacity to meet the total demands of served DCs and PCs. Constraints (17) and (18) indicate that the total route length of truck s to visit DCs or PCs should be less than or equal to planning horizon Ls. Constraints (19) and (20) ensure the total number of trucks serving DCs and PCs. Constraint (21) stipulates that each DC can be served by truck s within its time windows. Constraints (22) and (23) represent the departure and arrival of trucks at DCs, respectively. Constraint (24) stipulates that each PC can be served by truck s within its time windows. Constraints (25) and (26) represent the departure and arrival of trucks at PCs, respectively.
Constraints in the second echelon: 1, 1, , , , Constraint (27) indicates that each vehicle can depart from only one facility. Constraint (28) stipulates that each delivery customer i can be served by a single vehicle from DC or delivery customer m within the kth period. Constraints (29) and (30) ensure that each pickup customer j can be served by a single vehicle from delivery customer i or from PC and pickup customer n within the kth period. Constraints (31)- (33) imply that each vehicle leaves to serve another customer after serving one customer. Constraints (34) and (35) can be used to eliminate subtours on every route departing from each DC or PC to serve customers. Constraints (36) and (37) guarantee that each vehicle has sufficient delivery or pickup capacity to meet the total demand of served delivery or pickup customers. Constraint (38) states that the total route length in visiting the customers for vehicle v should be less than or equal to the planning horizon Lv. Constraint (39) ensures the total number of vehicles serving customers. Constraint (40) stipulates that each delivery customer i can be served by vehicle v within its time windows. Constraints (41) and (42) represent the departure and arrival of vehicles at delivery customer i, respectively. Constraint (43) states that each pickup customer j can be served by vehicle v within its time windows. Constraints (44) and (45) express the departure and arrival of vehicles at pickup customer j from delivery customer i, respectively.
Constraints (46) and (47) represent the departure and arrival of vehicles at pickup customer j from PC or pickup customer n, respectively.

Constraints in the first and second echelons:
, , Constraint (48) stipulates that the demand quantity of DC d is equal to the total demand quantities of assigned delivery customer i. Constraint (49) shows that the demand quantity of PC p is equal to the total demand quantities of assigned pickup customer j. Constraint (50) represents that vehicle v should arrive at DC d within its time windows. Constraints (51) and (52) define the departure and arrival of vehicles at DC d, respectively. Constraint (53) expresses that vehicle v should arrive at PC p within its time windows. Constraints (54) and (55) guarantee the departure and arrival of vehicles at PC p, respectively.

Solution Methodology
As an extended generation of the VRPPD [29,30], the MDPVRPPD can be considered an NP-hard problem that can be solved using exact approaches to achieve routing optimization with instances up to small-scale customers under reasonable computation times [11,48]. Heuristic or metaheuristic approaches are commonly developed as solution methodologies for large-scale instances. The solution approach for the VRPSPD involves the following steps. First, customers are clustered in accordance with the route constraints. Second, vehicles are assigned to clusters for trips. Finally, the traveling salesman problem is solved in each cluster of customers [45]. In this study, a two-stage algorithm combining a 3D k-means clustering and an IMOPSO algorithm is used to determine the open-closed routes for pickups and deliveries and find Pareto optimal solutions with maximum diversity and convergence. The hybrid optimization algorithm is shown in Figure 2. Here, k is set as the number of clusters, v(0) and p(0) are the initial velocity and position of each particle, respectively. run denotes the number of optimization runs, rmax is the maximum number of optimization runs, t is the t-th iteration, and it_max is the maximum number of iterations.

Three-Dimensional (3D) Clustering
The efficiency of clustering in reducing the complexity of logistics networks and improving computational efficiency has well been demonstrated when dealing with large-scale logistics network optimization [49,50]. K-means and multidimensional clusterings have been used to classify customers on the basis of multiple characteristics and obtain initial solutions for optimal vehicle routes [51,52]. The 3D K-means clustering approach adopted in this paper focuses on reassigning each customer to their new nearest depot in accordance with the geographic coordinates and time windows, and generates the initial population of optimized vehicle routes for the next multiobjective optimization approach.
In the present study, the clustering procedure can be stated in Algorithm 1.

Algorithm 1 3D K-means clustering
Input: Coordinates and time windows of depots and customers Output: Clusters composed of depots and their assigned customers 1. Select three one-dimensional vectors, namely, longitudes, latitudes, and time windows from depots and customers; 2. Combine the vectors to construct the initial 3D data elements; 3. Determine p service periods through the time vector classification; 4. For i =1:p 5. Set the number k of clusters on the basis of the number of depots; 6. Select the data from depots as the initial center points of each cluster; 7. Traverse all the data elements; 8. Calculate the distances from each element to the center points; 9. Assign each element to the center point with shortest distance; 10. Update the k clusters and calculate the new center points; 11. If customers need to be adjusted among clusters; 12.
Return to Step 7; 13. Else 14. End if 15. End for 16. Output the clustering results As shown in Algorithm 1, the clustering process can be divided into distinct components: First, extract the geographic coordinates (i.e., longitudes and latitudes) and time windows of depots and customers to construct the initial 3D data elements. Second, determine the p service periods on the basis of time window categories. Third, set the number k of clusters on the basis of the number of DCs that can be selected as the initial k center points. Fourth, traverse all data elements and calculate the distances from each element to the center points. Fifth, assign each element to its nearest center point and the elements belonging to the same center point that can be observed as one cluster. Then, update k clusters and check if customers must be adjusted among clusters. If yes, return to Step 7, otherwise, output the clustering results. Finally, output the clustering results composed of depots and their assigned customers. The departure of a vehicle from each pickup or delivery depot to its customers in the multiperiod logistics network can be initialized according to the clustering results and optimized through the IMOPSO methodology.

IMOPSP
The PSO algorithm is a population-based search technology inspired by the feeding behavior of birds and fish and their methods of information exchange; here, each particle is regarded as a feasible solution with a specific position and velocity [53,54]. A multiobjective PSO (MOPSO) algorithm was developed on the basis of the Pareto front to solve the multiobjective optimization problems [55,56], and its efficiency in solving the multiobjective routing optimization problems has well been demonstrated [57]. When dealing with such problems, the related issues should be analyzed as follows: (1) how to obtain nondominated solutions and find Pareto optimal solutions; (2) how to update the particles through the personal and social best selection mechanisms; and (3) how to avoid local optimal conditions and maintain the convergence and diversity of solutions [58][59][60]. In this study, an IMOPSO algorithm integrating the archive-based selection mechanism is developed to improve the diversity and convergence of optimal solutions.

Standard MOPSO
In standard MOPSO, each particle among the swarm represents a potential optimal solution and can be updated until the optimal solutions are found during multiobjective optimization. We suppose that the swarm size is N. Each particle i (i ϵ N) contains a position vector pi and a velocity vector vi that can be determined by the number d and value range of variables. Each particle can be expressed as follows: pi and vi can be updated using the following equations: , , , where pbest defines the best previous position of the particle, and gbest is the best position among all particles of the swarm. The inertia weight w guarantees the impact of a previous velocity on the current velocity. cp and cg denote the acceleration weights of particles toward pbest and gbest, respectively, and rp and rg represent random values uniformly distributed within [0,1]. The second term of Equation (57) reflects the individual cognition of particles toward pbest and the third term signifies the social cognition of particles towards gbest among the swarm. The updated process of biobjective optimization is illustrated in Figure 3.  In Figure 3a, all the particles move towards the optimal positions, separately. Figure 3b shows that the movement of each particle with respect to its current position, pbest and gbest. Different from a single-objective optimization problem, the particles in the biobjective optimization cannot move in a single direction in f1 or f2. Instead, they tend to move towards an optimal-line that is called the Pareto front [53,55], which is specified in the following section.

Pareto Front in the IMOPSO
In single-objective PSO, the particle with the best position can be determined as the leader in the swarm. In biobjective or multiobjective optimization, a series of Pareto optimal solutions can be stored in the Pareto fronts rather than one global optimum. Multiple Pareto front solutions are not dominated by any other solutions and can be selected as candidates of the best solutions [54,59]. One type of Pareto front surface is found in triobjective optimization. A schematic of the Pareto front solutions is shown in Figure 4.  In Figure 4a, the Pareto front surface is constructed from the Pareto solutions with three objectives, on the basis of the Pareto front line in the biobjective optimization. In Figure 4b, there exist three types of particles in the MOPSO, the particles under, on, and above the Pareto surface. Each particle has its distance to the Pareto front surface, and each particle can be updated to shorten the distance by optimizing its velocity and position. During the updating process, the gbest can be selected from the particles on the Pareto front surface, and the improved selection mechanism in this study is shown in the next section.

Archive-Based Selection Mechanism
In the updated procedure of each particle formulated in Equation (61), each particle can be adjusted toward the new optimal position based on pbest and gbest. The leader (gbest) with the best position in the swarm can be randomly selected from the Pareto front solutions to improve the diversity of the optimal solutions. An improved archiving strategy with nondominated sorting [61], a geographical-based system (grids) [53,59], and crowding distance [62,63] is utilized in this study to ensure the diversity and convergence of optimal solution selection from the archive. All potential solutions can be decomposed into nondominated (Pareto-dominance) or dominated solutions after nondominated sorting. Nondominated sorting can be rapidly carried out via the following steps: Step 1: Evaluate the objectives of each solution and calculate the two parameters of each particle as np and sp. np defines the number of solutions dominating solution p, and sp is the number of solutions dominated by p.
Step 2: Find all solutions in which the value of np equals to 0, save them in the first nondominated level called the Pareto optimal front, and set their Rank as 1.
Step 3: For each solution p in Rank 1, traverse each member q in its sp and find the value of nq minus one. If the value obtained is 0, save it in the second nondominated level and set its Rank as 2.
Step 4: Apply the above procedure for each solution in Rank 2 to obtain the third nondominated level and set its Rank as 3. Repeat the process until all nondominated levels are identified.
The external archive can be updated on the basis of nondominated sorting and is illustrated in Algorithm 2. We suppose that the archive size in the current iteration is Nt and that the maximum archive size is Na.
Evaluate the objective functions of the particle 3.
If particle i is dominated by the particles in A(t-1)

4.
Particle i cannot join the new external archive 5. Else 6.
Delete the worst particle and add particle i to A(t) 8.
Delete the particles dominated by i and add i to A(t) 10. End 11. End 12. End The density of particles in the Pareto front in the updated archive can be calculated through the grids and crowding distance methods shown in Figure 5, and particles in the grid with the fewest particles or within a low-crowding area can be selected in a large rate from all particles in the best Pareto front. In Figure 5a, the explored search space is divided into several grids and the particle in the grid with the fewest individuals is selected as the leader with the best position. The gbest selected in the grid has only one member, and other grids contain more than one member. In Figure 5b, the distances among the optimal solutions can be calculated as follows: First, sort the particles in accordance with their objective functions. Second, find the particles with the smallest and largest values (fm min and fm max ) of each objective function m and assign them an infinite value of crowding distance. Finally, calculate the distance distm(i) of solution i in the function m of two neighboring solutions (i−1 and i+1) by using Equation (59).
Then, the crowding distance CDi of solution i can be calculated as the sum of distance distm(i) in m objective functions by using Equation (60).
The pbest of each particle in the swarm in the current iteration pbest(t) can be updated on the basis of its former best position pbest(t−1) and new position pt and calculated using Equation (61). The overall framework of the IMOPSO is shown in Algorithm 3. For i=1:part_size 6.

Algorithm 3 IMOPSO
Select gbest from archive A(t) using grids or crowding distance 7.
End 11. Calculate the objective functions of new particles 12. Update archive A(t) 13. End

Minimum Costs Remaining Saving (MCRS) for Profit Allocation
The MCRS proposed by Heaney and Dickinson [64] is an effective profit allocation method to reasonably obtain the allocated cost for each participant in the collaboration [65,66]. Allocated costs in the MCRS can be calculated on the basis of the participants' feasible minimum and maximum costs. We suppose that the set of all participants is T= {1, 2, …, t}, wtmin is the lower profit of participant t, and wtmax is the upper profit of participant t. L represents the subset of T. The profits of alliances T and L are determined as W(T) and W(L), respectively. W(L) can be calculated in Equation (62) as follows.
where ξ is the synergy coefficient of the organizer to coordinate the collaboration. C0(t) defines the original operational cost of participant t. C(L) denotes the total cost of participants in alliance L. φt (T,w) defines the allocated profit of participant t in alliance T and can be calculated in Equations (63) and (64) as follows.
where wtmin and wtmax can be obtained using the following equations: where w0(t) shown in Equation (65) is the profit of participant t from the initial network to the optimized network without joining the collaboration. W(L−{t}) shown in Equation (66) is the remaining profit of alliance L excluding participant t. Equation (67) guarantees that the total profits of participants are higher than the profit in alliance L. Equation (68) defines the lower bound of wtmin, while Equation (69) ensures the lower and upper bounds of φt (T,w). The profits of each collaborative member can be fairly allocated and the collaboration can be promoted.

Data Source
The effectiveness and applicability of the proposed multiobjective optimization model and hybrid optimization algorithm are verified through a real logistics network with pickup and delivery services in Chongqing City, which is considered one of the cities with the fastest development of the e-commerce platform and logistics transportation in China. The logistics network shown in Figure 6, includes three DCs (DC1, DC2, and DC3), three PCs (PC1, PC2, and PC3), and 175 customers (C1, C2, …, C175) for pickups and deliveries among the three service periods. The characteristics of the logistics facilities and customers are provided in Table 4.

Optimization Results
The MDPLNPD can be simplified through 3D k-means clustering described in Section 5.1 on the basis of the characteristics of facilities and customers shown above. Multiple service periods can be obtained, and each customer can be reassigned to an adjacent facility in accordance with the geographic coordinates and time windows. The clustering results are illustrated in Figure 7; here, the large spheres represent DCs and PCs and the small spheres denote customers. A large sphere and its neighboring small spheres of the same color can be set in the same cluster. In Figure 7, facilities and the customers with similar characteristics are grouped into one cluster on the basis of their geographic coordinates (i.e., longitudes and latitudes) and the center point of their time windows. Three service periods are determined through the three categories of time windows. Then, the MDPLNPD can be reconstructed and the travel routes from facilities to customers in each service period can be optimized using the proposed IMOPSO algorithm in Section 5.2. We suppose the following parameters for the IMOPSO algorithm: w = 0.9, cp = 1, rp = 2, part_size = 200, It_max = 150, and Na=100. The optimization results of total operational cost (Cost), waiting time (Time), and number of vehicles (Vehicle) in the 1st-2nd-3rd period collaboration are shown in Table  5. As shown in Table 5, the gaps in total operational cost, waiting time, and number of vehicles indicate cost savings from the initial network to the optimized network obtained through collaboration and resource sharing. For instance, the cost of the alliance {PC1, PC2, PC3} decreases by 21.5%, the time reduction of {DC1, DC2, PC1, PC2} is 48 min, and the number of vehicles required by the alliance {DC1, DC2, DC3, PC1, PC2} decreases to 16. The reduction percentages of gaps increase when collaboration exists because of the addition of new members. For example, when PC3 joins the existing alliance {DC1, DC3, PC1, PC2}, the gaps of operational cost, waiting time, and number of vehicles increase by 38.8%, 23.5%, and 23.1%, respectively. Moreover, the grand alliance {DC1, DC2, DC3, PC1, PC2, PC3} achieves the maximum gaps of operational cost, waiting time, and number of vehicles as $2824, 75 min, and 19, respectively. Therefore, the grand alliance {DC1, DC2, DC3, PC1, PC2, PC3} can be selected as the best alliance in the collaborative MDPLNPD. The optimized vehicle routes from collaborative facilities DC1 and PC1 in the grand alliance are described in Table 6 and Figure 8; in this collaboration, the vehicles are shared to reduce their number and maximize resource utilization.   In Table 6 and Figure 8, the vehicles can be shared among facilities at different service periods to improve transportation efficiency and optimize the resource configuration. In Figure 8a, the first vehicle departs from DC1 to visit delivery customers, visit pickup customers, and return to PC1 in the first period. In Figure 8b, the second vehicle used for pickups by PC1 in the first period can be shared by DC1 to achieve deliveries in the third period. In Figure 8c, the open route (e.g., DC1→C96 →C110→C142) for deliveries by DC1 in the second period and the open route (e.g., C40→C28→C27 →C23→PC1) for pickups by PC1 in the third period can be operated by a third vehicle because of vehicle sharing. Therefore, the resource can be reasonably allocated on the basis of several types of vehicle sharing to promote effective collaboration among multiple facilities and service periods. A comparison of the results of pickup cost (PC), delivery cost (DC), maintenance cost (MC), penalty cost (Pen_cost), total cost (TC), waiting time (WT), and number of vehicles (Vehicle) from the noncollaborative network to the collaborative network in the grand alliance is shown in Table 7.  Table 7 shows that the total operational cost, waiting time, and number of vehicles decrease on the basis of collaboration and resource sharing in the MDPLNPD. For example, the pickup and delivery costs are reduced by 22.2% and 25.6%, respectively. The number of vehicles decreases from 30 to 11 because of the vehicle sharing among multiple facilities and periods; thus, the maintenance cost also decreases from $450 to $154. The waiting time is reduced from 111 min to 36 min, and the penalty cost reduces from $555 to $180. Consequently, effective collaboration and resource sharing among multiple facilities and periods can be used to enhance transportation efficiency and maximize resource utilization during open-closed mixed network optimization for the MDPLNPD.

Analysis and Discussion
We divide the entire study period into three collaborative periods. Facilities can collaborate within each collaborative period or across the three collaborative periods. We, therefore, consider four scenarios: in each of the first three scenarios, facilities collaborate only within each of the three collaborative periods; and in the fourth scenario, facilities also collaborate across the three collaborative periods. Multiobjective optimization in the MDPLNPD is performed for these four scenarios and the calculation results are shown in Table 8 and Figure 9.  1st  2689  16  5  3296  37  10  2nd  2331  14  4  3162  34  8  3rd  2804  17  6  3589  40  12  Total  7824  47  15  10047  111  30  1st-2nd-3rd  7223  36 11 . Figure 9. Illustration of optimization result comparison among the four scenarios. Table 8 and Figure 9 demonstrate that the total operational cost, waiting time, and number of vehicles in the 1st-2nd-3rd period collaboration achieve the lowest values compared with the total in the other scenarios of collaboration. The operational cost decreases by 18.4%, 26.2%, 21.8%, and 28.1% in the four scenarios from the initial network to the optimized network. The first, second, third, and 1st-2nd-3rd period collaborations achieve reductions in the number of vehicles of 5, 4, 6, and 19, respectively. The waiting time of the 1st-2nd-3rd period collaboration has the minimum value of 36 min, compared with the total waiting time of 47 min in the first-, second-, and third-period collaborations. Therefore, the MDPLNPD optimization (across the three collaborative periods) outperforms the collaboration in each of the three collaborative periods. The collaborative stability can be ensured on the basis of the profit allocation method in Section 5.3. The profit allocation in the grand alliance {DC1, DC2, DC3, PC1, PC2, PC3} in the four scenarios can be calculated, and the allocated cost savings of each collaborative member are shown in Table 9 and Figure 10.   Table 9 and Figure 10 demonstrate that members in the 1st-2nd-3rd period collaboration achieve the largest allocated cost savings {475, 500, 513, 454, 449, 433} than in the other three scenarios. The total allocated cost savings in the first scenario (first-period collaboration) have the minimum values {95, 83, 106, 125, 111, 87} for each participant. The total cost savings for each participant in the first-, second-, and third-period collaborations are lower than that in the 1st-2nd-3rd period collaboration, which indicates that each facility would likely choose to join the 1st-2nd-3rd period collaboration for maximum cost savings. Therefore, the grand alliance formed with the 1st-2nd-3rd period collaboration can be selected as the best collaborative strategy for the logistics enterprises, which provides the highest benefits and strongest stability.

Management Insights
In this study, a collaborative mechanism is effectively used to coordinate the open and closed mixed vehicle routes in a multidepot logistics network with pickups and deliveries. Resource sharing strategies among multiple service periods are introduced to promote the sustainable development of city logistics. The management insights gained from this paper are summarized as follows: 1. Collaboration and resource sharing strategies are developed to optimize the MDPVRPPD and achieve improved the operational efficiency and resource configuration. A collaborative mechanism is designed using multiple facilities and periods to coordinate open and closed routes with enhanced operational efficiency; here, pickup and delivery services can be incorporated to expand the practical applications of the system. Resource sharing strategies combining customer information, transportation resources, and facility services are investigated in the MDPLNPD to achieve maximum resource utilization and maintain collaborative stability. The results of the proposed mathematical model and solution methodology serve as a reliable reference for the applicability and validity of our collaborative mechanism and resource sharing strategies in addressing the MDPVRPPD. Collaboration and resource sharing are important factors in ensuring a reasonable resource configuration to achieve the sustainable development of smart city logistics. 2. The collaborative mechanism and resource sharing strategies applied for multidepot multiperiod logistics network optimization can guide management decisions for enterprises to survive in the competitive logistics market and reduce their environmental impact to achieve sustainable logistics systems. Considering that sustainability is an important societal concern requiring effective and flexible management, the design and operation of a sustainable logistics system are then a challenge faced by many enterprises. The adoption of new information technologies into logistics systems can produce a reasonable logistics resource configuration that can ensure an intelligent logistics system. Long-term collaborations among multiple logistics enterprises should be maintained and resource sharing strategies should be promoted by the concerned governments to construct a resource-friendly logistics system and propel sustainable logistics development.

Conclusions
Collaboration and resource sharing are introduced in this study to optimize the MDPLNPD. Specifically, a collaborative mechanism is established to incorporate the open-closed mixed vehicle routes with enhanced operational efficiency and the resource sharing strategies are investigated to optimize the resource configuration with the maximum resource utilization. This study has three major contributions. First, a multiobjective optimization model is developed to formulate the collaborative MDPLNPD with minimum operational cost, waiting time, and number of vehicles. Next, a hybrid heuristic algorithm combining 3D k-means clustering and the IMOPSO algorithm is proposed to solve the MDPVRPPD. Here, 3D k-means clustering is performed to simplify the network structure and reduce the waiting time. The IMOPSO algorithm is employed to find Pareto optimal results. Third, a profit allocation method is developed to distribute the cost savings and ensure the stability of the collaboration of each collaborative member.
A case study in Chongqing City, China is analyzed to demonstrate the performance of the proposed multiobjective model and solution methodology in improving collaboration and resource sharing in MDPVRPPD. The numerical results show that the total operational cost, waiting time, and number of vehicles could be reduced by $7223, 36 min, and 11 in the optimized network, respectively. A comparison of the performance of several algorithms shows that the IMOPSO used in this study is superior to other algorithms, such as NSGA-ΙΙ and MOEA for finding high-quality Pareto optimal solutions. The MCRS method is used to allocate the cost savings to each participant within a collaborative alliance, which maintains the collaborative stability. Furthermore, comparison among the four types of period collaborations shows that the period collaboration involving all three time periods achieves the minimum operational cost, waiting time, and number of vehicles. The allocated cost saving of each collaborative facility in the optimal period collaboration type is significantly higher than other types of period collaborations. We thus identify the best collaborative strategy with the three service periods in the case study.
From a practical perspective, this research encourages collaboration and resource sharing to improve the robustness and reliability of the collaborative MDPLNPD. On the one hand, resource sharing strategies are proposed to improve the utilization of transportation resources and the efficiency of facility service. This strategy can lead to an effective resource configuration and serve as a decision tool for urban logistics departments. On the other hand, the combination of the collaborative mechanism and resource sharing strategies leads to the formulation of a grand alliance in the collaborative MDPLNPD, which proves more beneficial than the noncollaborative MDPLNPD and the scenarios of single period collaboration. In addition, customer demand characteristics and the heterogeneity of collaborative alliances will affect the stability and reliability of the collaborative MDPLNPD. Therefore, the collaborative mechanism and resource sharing strategies can be further studied and used as a reference for logistics companies and local transportation departments to establish sustainable intelligent logistics systems.
This study introduces applicable and operational methodologies to solve the MDPVRPPD with collaboration and resource sharing, and provides logistics managers with practical strategies. The limitations of the present work and further research along the line in the future can proceed in the following directions. (1) The influence of demand fluctuations and customer satisfaction can be considered in the MDPVRPPD and studied how to achieve collaboration in the MDPLNPD. (2) Dynamic mathematical models and exact algorithms can be adopted and integrated into hybrid heuristic algorithms to determine more accurate solutions for the MDPLNPD. Acknowledgments: The authors would like to express our sincere appreciation for the valuable comments made by three anonymous reviewers, which helped us to improve the quality of this paper.

Conflicts of Interest:
The authors declare no conflicts of interest.