Network Planning Method for Capacitated Metro-Based Underground Logistics System

Abstract

Underground logistics system (ULS) tends to alleviate traffic congestion, increase city logistics efficiency, mitigate the negative effects of traditional logistics processes, and improve the sustainability of urban areas. However, the relatively high cost and risk of underground construction are serious obstacles to implementing ULS. Integrating ULS into modern metro system (M-ULS) is considered to be feasible and efficient to solve this problem. This paper aims at developing a metro system-based ULS network planning method. First, an evaluation model of underground freight volume was proposed considering service capacity, freight flow, and regional accessibility. Second, a set of mixed integer programming model was developed to solve the problem of optimal nodes’ location-allocation (LAP) in the network. Then, a hybrid algorithm was designed with a combination of E-TOPSIS, exact algorithm, and heuristic algorithm. Finally, two lines of Nanjing Metro were selected as a case to validate the proposed planning method. The results showed that the new system can significantly reduce the construction costs of ULS and alleviate traffic congestion. Moreover, the potential of metro stations and underground tunnels can be fully exploited to achieve higher logistics benefits.

This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

1. Introduction

Traffic congestion is one of the “dilemmas” encountered by metropolises around the world. The fifth type of transportation and supply system, “underground logistics system” (ULS), has been developed to address the energy consumption, environmental pollution, and safety issues [1]. Various forms of ULS with different locomotives or driving methods have been developed and applied, such as CargoCap underground freight transportation system in Germany [2], OLS-ASH project in Netherland [3], and Pipe§net vacuum freight capsule in Italy [4–7]. ULS is widely recognized as a great innovation in the logistics industry because of its huge benefits to the city. It can quickly transport goods through underground pipes or tunnels within cities or between cities [8, 9]. However, it is undeniable that the high cost, long construction period, and high-risk of underground projects [10–12] are serious obstacles to ULS development. Reducing the construction cost is critical to improve the uptake of ULS [13].

Dampier and Marinov proposed a new conception of “integrating underground freight transport into modern metro system” (metro-based ULS, designated as M-ULS), which offered a feasible way with limited budget [14]. As a special form and promising alternative of ULS network, M-ULS has become one of the hottest research topics. Especially, setting up the intermodal network aims at maximizing the efficiency of cargo transport without interfering with the passenger transport system. The organizational process of this new system is complicated. It involves the linkage and optimization of the dual objectives of the passenger and cargo transportation system. The cargo needs to be sent to the demand points via reasonable metro stations and the distribution system on the ground. Therefore, the global optimization of the M-ULS efficiency can be classified as a special type of capacitated multistage facility location problem (CMFLP). The aim of this paper is thus to design a multilevel node location-allocation solution for M-ULS through mathematical modeling and simulation.

Research on ULS has been ongoing for nearly 30 years, but M-ULS is still an emerging area. Existing studies verified the feasibility of M-ULS and proposed concept designs. M-ULS is defined as an intermodal mechanism of mixed passengers and freight, in which batch driving or attachment carriages are enabled by utilizing metro railways [15]. Their advantage not only lies on the reservation of traditional ULS benefits in alleviating traffic congestion, improving urban environment and freight efficiency, but also contributes to a huge saving of the underground space, reducing construction costs and period. The passenger transport efficiency would not be interfered through rational technology design and operation management [15]. In Tokyo, new subway-integrated ULS received public endorsement in a pilot project [16]. In Newcastle upon Tyne, metro network was used for the collection of small-sized to medium-sized parcels, low-density high-value goods, and recyclable material [17]. In Korea, the freight volume of Seoul Metro 50 platforms was predicted to validate the reduction of social cost and environmental problem [18]. Toronto subway has also been evaluated on how such a modal shift appealed to cargo shippers [19]. All in all, real-world simulations were widely developed in recent years [20–22].

However, few articles with the topics on the planning method of M-ULS network have been published. Hörl et al. and Masson et al. conducted qualitative analyses of node location problem of single-line metro freight transport [23, 24]. Fangting et al. proposed a location-routing model of multiline transfer problem of M-ULS [25], which aimed at minimizing delivery time to improve the efficiency of metro freight. The optimization of the single-layer network topology was carried out; however, the influence of urban freight flow distribution on the system service level was not involved. Ambrosino and Scutella [26] proposed a four-level complex network distribution structure including facility location, warehousing, transshipment, and strategic decisions. The findings have contribution for analyzing the CMFLP of M-ULS, although it does not directly discuss M-ULS. M-ULS itself is a complex system. Functional classification of different types of facilities is necessary to meet the needs of intermodal transport. In order to pursue higher logistics benefits, a systematic evaluation of the urban freight transport demand is required due to strict capacity constraints of M-ULS. The freight flow direction and node location are two key indicators of M-ULS network layout. Thereby, this paper proposed a hierarchical mixed integer programming model (MIPM) combined with comprehensive evaluation that could reflect the different relationships of nodes. Compared to the modeling benchmark of general facility location problem [27–29], three additional important factors had been considered: (i) Constraint of metro network. Compared to the traditional logistics LAP model, M-ULS distribution path, which is restricted by the existing metro network, cannot be freely selected and adjusted immediately. (ii) New system service area. The definition of new system service objects has to be integrated into the optimization process because maximizing the limited system service capability is a key to M-ULS planning. (iii) Mixed transport organization is complicated. The design principle of existing metro network is passenger oriented, and the passenger transportation efficiency cannot be influenced by the new system. Improving the service performance of the new system requires redesigning an organizational process of passenger and freight mixed transportation, and the global optimization of freight transportation process.

The contributions of this study mainly lie in three aspects: (1) an effective organizational planning method for M-ULS was proposed; (2) a mixed integer programming model with multiple constraints was developed to solve the capacitated metro-based underground logistics system location-allocation problem (CM-ULSLAP); (3) a hybrid heuristic algorithm was developed for the simulation of real-world. Two lines of Nanjing Metro were selected as a case to validate the proposed planning method. The findings can promote the implementation of ULS. In addition, this study offers new ideas and quantitative optimization approaches for linear or network-based complex engineering projects. The supplied case study provides additional insights for research and practices.

2. Problem Description

The CM-ULSLAP can be described as follows: there are m operating metro lines, n urban distribution centers (designated as UDCs), and k demand points. Firstly, cargo from the open UDCs which provide ULS services is imported into underground network through the nearby metro stations (designated as ULS terminal). Then, dedicated cargo trains travel between ULS terminals and different underground depots (designated as UDs). Those UDs refer to the retrofitted metro stations, in which UD-oriented trains are bypassed to the parallel freight platforms before entering the passenger platforms. In this process, the frequency of freight fleet keeps interlaced with metro, and existing transfer points can be activated as underground hubs (designated as UHs) to support underground transshipment. Finally, freight units floating up from UDs are allocated to the corresponding ground terminals (designated as GTs) by ground electric vehicles. Cargo is delivered to the customers who are within the service scope of M-ULS after automated sorting and packaging at GTs. The reversed procedures are presented in mirror image. Figure 1 demonstrates the organizational breakdown structures with planning issues of M-ULS.

[figure omitted; refer to PDF]

The initial data that can be used for the metro freight planning usually focus on several incomplete customers characteristics, such as location, regular delivery route, and handover period record. Other important input data are the demand for travel, which is commonly formulated as the origin-destination (OD) matrix [30]. The freight OD matrix in this study is extracted based on the counts of travels multiplied by order capacity from one area to another. On this basis, the alternative solution aims to improve the current customer service with low cost and high efficiency. The high efficiency of M-ULS depends on the coordination of its own technical integration and freight flow distribution. In order to get this, the joint freight function should not lead to any degradation or disruption of passenger service [15]. As a result, goods and passengers must be separated. Meanwhile, the efficiency loss caused by intermodal transport can also be compensated by a highly integrated urban logistics benchmark. Taking into account the regional development density and spatial impact on the surrounding environment, in M-ULS, freight flows cannot be detained too much after reaching corresponding underground destinations. Instead, they should be transferred directly to the newly built ground terminals for processing. By setting such a transit layer between customers and underground distribution network, the transport efficiency and system service capacity can be largely improved by dividing the M-ULS from the traditional logistics systems. Based on the planning principles discussed above, the node facilities are divided into the following four layers:

2.1. ULS Terminal

The ULS terminal plays a role in docking UDC with underground network. The freight units which have been well sorted and packed in UDC will be simply coded at the ULS terminal to each exact ground terminal or underground depot. In this paper, the location of the ULS terminal was determined by calculating the weighted distance to nearby urban distribution centers.

2.2. Underground Depot (UD)

The UD is defined as an enhanced metro station which has the function of cargo handling and connects the underground logistics network with upper ground terminals. The number and distribution of UD directly affect the system service. In this paper, a MILP method with corresponding heuristic algorithm was introduced to search the optimal location-allocation of UD-GT.

2.3. Ground Terminal (GT)

GT is one of the optimization targets, which is affiliated with UD. GT is also a key part of the entire M-ULS supply chain because it is directly oriented to customers whose delivery requirements are time efficient. The cargo handling capacity of GT is the major handicap in the performance of system service. In this paper, the location of ground terminal was transformed into the capacitated set covering problem (CSCP). And a set coverage model for GT location is established.

2.4. Demand Point

Demand point is also one of the optimization targets. Because of the limited transportation capacity of the metro system which is far from meeting cargo supply for all customers, a screening evaluation of demand points is essential for maximizing logistics benefits and service scope.

In order to simplify the study, the following assumptions for CM-ULSLAP are proposed: (1) suburban UDC reasonably allocates goods into the underground network, and the freight handling capacity of ULS terminal is unlimited; (2) each candidate GT can only correspond to one UD, and the source of freight flows is unique; (3) taking a full carriage of dedicated trains as a freight unit, the less than truckload (LTL) transportation are not allowed in the underground distribution network; (4) the processing time and transfer cost at underground hub are not considered, and all UHs are activated by default; (5) intermodal transport does not involve the urban road network layout and additional expenses caused by traffic congestion; (6) M-ULS does not reschedule the original timetable of metro system, but underground freight flow frequency can be adjusted based on the traffic volume of passengers; (7) construction and storage condition at UD are not considered; and (8) only part of customers’ attributes are available.

3. E-TOPSIS for Evaluation of Underground Freight Volume

The entropy weighted TOPSIS method (E-TOPSIS) is applied in the weighting process of the evaluation indicator system [31]. The information entropy calculated in the indicator matrix is used to weight the dispersion degree of demand characteristics, all that can avoid the influence of human subjective factors [32]. The final points of the evaluation objects in multiobjective decision-making are gotten by the TOPSIS method which can solve the higher requirement for the sample.

This paper used unity of freight flows (UFFs), regional accessibility (RA), and source distance (SD) as three dimensions to reflect UDC’s preference for the characteristics of customers’ demand. Therein, variables UFF and RA were obtained from Nathanail et al. [33], who maintained that the urban traffic alleviation and the development of integrated transport to reduce the scattered freight delivery were two of the most remarkable characteristics of the sustainable city logistics. In our study, we took these two quantitative features as the functional requirements of the metro freight system, where UFF was formulated as the accumulated OD from different directions according to Masson et al. [24] and RA was formulated as the accumulated delivery period from different directions divided by average path length according to Montoya-Torres et al. [34]. SD was formulated as a variable that supplemented UFF with distance influence, which involved the data of spatial distance and freight OD from different directions. These 3 variables share the initial data mentioned above.

(i) Unity of freight flows. Since UD does not have the function of parcel disassembly, logistics benefits are promoted along with the growth of underground mileages [35]. Therefore, M-ULS should provide centralized supplies as much as possible to areas with huge demand:$$\begin{array}{}\text{(1)}& x_{1j}=\overline{Q_j},\end{array}$$where $Q$ is the total amount of freight picked up and delivered per day at demand point $j$.

(ii) Regional accessibility. The efficiency of urban logistics is largely influenced by traffic congestion [36–38]. By removing goods from above to underground, the freight mobility of some regions with poor accessibility can be effectively improved [39]. The index of regional accessibility can be characterized by the time-road response function:$$\begin{array}{}\text{(2)}& x_{2j}={\displaystyle \sum }_{i=1}^m{\left(\frac{t_{ij}}{S}\right)}^{{\tau }_{ij}},\end{array}$$where $t_{ij}$ is the average delivery time form UDC i to demand point j, $S$ is the length of the delivery route, and ${\tau }_{ij}\in \left(\mathrm{0,1}\right)$ is the time sensitive factor.

(iii) Source distance. Underground logistics can generate considerable internal and external benefits [40]. M-ULS should take full advantage of the length of metro lines. Putting goods which are far from the source (ULS terminal) but in huge demand into the system will lead to the overall increase of logistics values. The source distance index can be portrayed as$$\begin{array}{}\text{(3)}& x_{3j}={\displaystyle \sum }_{i=1}^m‖d_{ij}‖\cdot {\left[{\gamma }_i\left(x_j\right)\right]}^{{\alpha }_j},\end{array}$$where $‖d_{ij}‖$ is the Euclidean distance between UDC i and demand point j, ${\gamma }_i\left(x_j\right)$ is the freight correction factor of the sample, which represents the proportion of the freight volume from direction i, and ${\alpha }_j$ is an amplification coefficient.

An evaluation model for the priority level of M-ULS service is constructed.

(1) The initial decision matrix. Assuming there are m evaluation indicators and evaluation indicator sets have n subsets, the evaluation value of indicator i in subset j is $x_{ij}$. The decision matrix of all subsets is$$\begin{array}{}\text{(4)}& R={\left(x_{ij}\right)}_{m\times n}=\left[\begin{array}{ccc}{x}_{11}& \cdots & x_{1n}\\ ⋮& \ddots & ⋮\\ x_{m1}& \cdots & x_{mn}\end{array}\right].\end{array}$$

(2) Standardization of the decision matrix. In order to solve the uniform of indicators’ units, this paper makes the standardization on all indicators. Take the best value of indicator i as $x_i^{\text{best}}$, the worst as $x_i^{\text{worst}}$. The efficacy coefficient $\alpha$ is introduced to obtain the normalized data matrix B:$$\begin{array}{}\text{(5)}& B={\left[x_{ij}\right]}^{\prime }=\frac{x_{ij}-x_i^{\text{worst}}}{x_i^{\text{best}}-x_i^{\text{worst}}}\cdot \alpha +\left(1-\alpha \right),\hspace{1em}\alpha \in \left(\mathrm{0,1}\right).\end{array}$$

(3) Determine the indicators weight. The indicator weights are determined by the expert evaluation method in the TOPSIS. The results are quite subjective. Therefore, this paper adopts the information entropy method to determine the indicators weight.

$e_i$ denotes the entropy of indicator i in the standardized decision matrix:$$\begin{array}{c}\text{(6)}& T_{ij}={\displaystyle \sum }_{k=1}^l\left(\frac{x_{ij}}{{{\displaystyle \sum }}_{j=1}^n{x}_{ijk}}\right),e_i=h\cdot {\displaystyle \sum }_{j=1}^n{T}_{ij}\ln \left(T_{ij}\right),\hspace{1em}\forall i=\mathrm{1,2},\dots ,m,\\ h=\frac{1}{\ln n},\hspace{1em}0\le e_i\le 1.\hspace{1em}\text{If}\hspace{1em}{T}_{ij}=0,\hspace{1em}{T}_{ij}\ln \left(T_{ij}\right)=0.\end{array}$$

The dispersity of the evaluation value of indicator i can be represented as follows:$$\begin{array}{}\text{(7)}& {\beta }_i=1-e_i,\hspace{1em}\forall i=\mathrm{1,2},\dots ,m.\end{array}$$

If $T_{ij}$ is more dispersed, ${\beta }_i$ is bigger and indicator i is more important. If $T_{ij}$ is relatively concentrated, indicator i is less important. If all $T_{ij}$ are equal and the distribution is absolutely concentrated, indicator i is invalid. So, the weight factor of nonsubjective preference for indicator i is$$\begin{array}{c}\text{(8)}& Y={\left(y_1,y_2,\dots ,y_m\right)}^T,\forall y_i=\frac{{\beta }_i}{{{\displaystyle \sum }}_{i=1}^m{\beta }_i}.\end{array}$$

(4) The weighted matrix of indicators’ value. According to the normalized data matrix B and the obtained entropy weight $y_i$, the weighted matrix of indicators’ value is calculated by following formula:$$\begin{array}{}\text{(9)}& K=Y\cdot \left[\begin{array}{ccc}{T}_{11}& \cdots & T_{1n}\\ ⋮& \ddots & ⋮\\ T_{m1}& \cdots & T_{mn}\end{array}\right]=\left[\begin{array}{ccc}{k}_{11}& \cdots & k_{1n}\\ ⋮& \ddots & ⋮\\ k_{m1}& \cdots & k_{mn}\end{array}\right].\end{array}$$

(5) Calculating the relative distance. The selection scheme of candidate GT for this phase is calculated with the relative proximity ${\rho }_j$ of positive and negative ideal solutions under weighted matrix K:$$\begin{array}{}\text{(10)}& d_j^{+}={\left[{\displaystyle \sum }_{i=1}^m{y}_i{{\displaystyle \sum }_{k=1}^l\left(k_{ij}-x_i^{+}\right)}^2\right]}^{1/2},\hspace{1em}\forall j=\mathrm{1,2},\dots ,n.\end{array}$$where $d_j^{-}$ denotes the Euclidean distance between $k_{ij}$ and negative ideal solution:$$\begin{array}{c}\text{(11)}& d_j^{-}={\left[{\displaystyle \sum }_{i=1}^m{y}_i{{\displaystyle \sum }_{k=1}^l\left(k_{ij}-x_i^{-}\right)}^2\right]}^{1/2},\hspace{1em}\forall j=\mathrm{1,2},\dots ,n,{\rho }_j=\frac{d_j^{-}}{\left(d_j^{+}+d_j^{-}\right)}.\end{array}$$

For the candidate location scheme, ${\rho }_j$ is generally between 0 and 1. The closer to 1, the more appropriate the corresponding scheme. Finally, demand points in the newly built GT are selected to accept the M-ULS service according to the value of ${\rho }_j$ listing in descending order, until one of the nodes or underground sections is overloaded. Underground freight OD can be uniformly formulated as$$\begin{array}{c}\text{(12)}& \overline{Q}=Q_i^x\cdot A_j,\forall A_j=\left\{\begin{array}{c}{\left(\mathrm{1,1},\dots ,1\right)}^T,\hspace{1em}1\ge {\rho }_j\ge {\overline{\rho }}_j,\\ {\left(\mathrm{0,0},\dots ,0\right)}^T,\hspace{1em}0\le {\rho }_j\le {\overline{\rho }}_j.\end{array}\right.\end{array}$$

4. Multistage Node Location Model for M-ULS

4.1. Parameter Definition

(1)

Define $V_{\mathrm{C}},V_{\mathrm{S}},\text{\hspace{0.17em}and\hspace{0.17em}}{V}_{\mathrm{P}}$ as sets of metro station, candidate UD, and GT, respectively.

The cost parameters: $c_1$, integrated unit price for underground transport; $c_2=\beta \cdot c_1$, comprehensive price of ground transport, $\beta$ is the freight price coefficient, which reflects the pricing gap between two modes; and $c_3$, construction cost of UD (depreciated to daily).

The distance parameters: $d_{ij}$, the distance between UDC j and metro station i; $r_{ik}$, the shortest path for goods routed to the UD k, in which underground transfer is included; and $Z_x,\text{and}{Z}_i$, coordinates of GT and metro station.

The freight parameters: $Q_{ij}$, the amount of goods sent by UDC j to certain metro line where the UD i is located; $Q_k$, the amount of goods picked up by the UD k per day; $u_{\max }^j$, capacity of metro line j.

(2)

Decision variables.

$M_i$ : 1, if metro station i is selected as UD; 0, otherwise.

$P_{xk}$ : 1, if GT x is visited by the vehicles from UD k; 0, otherwise.

$H_{ij}$ : 1, if metro station i is selected as underground access of goods from UDC j; 0, otherwise.

4.2. Set Coverage Model for Ground Terminal

Sets $A\left(j\right)$ of all demand points covered by ground terminal $a_j$ and sets $B\left(i\right)$ of all GTs that cover demand point $b_i$ are given. Let $Y_{ij}$ be the distribution coefficient of cargo throughput for $a_j$ to $b_i$ and $X_j$ be the decision of whether to locate ground terminal at $a_j$:$$\begin{array}{}\text{(13)}& \min {\displaystyle \sum }_{j,a_j\in V_R^{\ast }}{X}_j.\end{array}$$subject to$$\begin{array}{}\text{(14)}& L\left(i,j\right)\le d,\hspace{1em}\forall b_i\in A\left(j\right);a_j\in B\left(i\right),\text{(15)}& {\displaystyle \sum }_{i,b_i\in A\left(j\right)}\overline{Q_i}{Y}_{ij}\le Q_j^{\max }{X}_j,\text{(16)}& {\displaystyle \sum }_{j,a_j\in B\left(i\right)}{Y}_{ij}=1,\hspace{1em}\forall b_i\in V_R^{\ast },\text{(17)}& X_j\in \left\{\mathrm{0,1}\right\},\hspace{1em}\forall a_j\in B\left(i\right).\end{array}$$

The objective (14) is to search the minimum number of GT under the current demand scenario; constraints (15) and (16) indicate the maximum service scope and cargo handling capacity of the ground terminal; and constraint (17) specifies the sum of freight demand $b_i$ shared by each GT should be equal to the amount of OD obtained in (13).

For unified presentation, the above model is written in the form of matrix inequality:$$\begin{array}{}\text{(18)}& \end{array}$$

K is a $2n\ast n\left(n+1\right)$ dimensional matrix. In the first n line, elements of column 1 to n equal 0; n + 1 to 2n is a diagonal matrix with principal diagonal $D_{11},\dots ,D_{n1}$, and so forth. In the line n + 1, column 1 to n principal diagonal is $-Q_{\max }^j$, while column n + 1 to $2n^2$ has a single line element $D_{11}\overline{Q_1},\dots ,D_{n1}\overline{Q_n}$, and so forth.

Z is a $n+n^2$ dimensional vector. b is 2n in length, of which elements in the first n line equal −1, otherwise 0. All solutions of Z satisfying $K\cdot Z\le b$, $0\le Y_{ij}\le 1$, $i,j\in V_R^{\ast }$ at the same time compose the feasible space of the GT location under the maximum service capability.

4.3. Location-Allocation Model for Underground Depot

$$\begin{array}{c}\text{(19)}& \min \hspace{1em}{J}_b={\displaystyle \sum }_{i\in V_C}{\displaystyle \sum }_{k\in V_S}{\displaystyle \sum }_{x\in V_P}{M}_i{P}_{xk}{Q}_x‖Z_x-Z_i‖c_1+{\displaystyle \sum }_{i\in V_C}{M}_i{c}_3+{\displaystyle \sum }_{i\in V_C}{\displaystyle \sum }_{j=1}^m{\displaystyle \sum }_{k\in V_S}{H}_{ij}\left(d_{ij}{Q}_{ij}\beta +r_{ik}{Q}_k\right)c_1,\end{array}$$ subject to$$\begin{array}{}\text{(20)}& Q_k\le Q_{\max }^k,\hspace{1em}\forall k\in V_S\text{(21)}& {\displaystyle \sum }_{j=1}^m{\displaystyle \sum }_{i\in V_C}{H}_{ij}{Q}_{ij}+{\displaystyle \sum }_{x=1}^q\left({\displaystyle \sum }_{k=x+1}^{y-1}{H}_{ij}{Q}_{xk}\right)\le u_{\max },\hspace{1em}\forall q\in \left\{\mathrm{0,1},\dots ,m-1\right\},\text{(22)}& {\displaystyle \sum }_{i\in V_C^n}{H}_{ij}=1,{\displaystyle \sum }_{j=1}^m{H}_{ij}\ge 1,\text{(23)}& {\displaystyle \sum }_{k\in V_S}{P}_{xk}=1,{\displaystyle \sum }_{x\in V_P}{P}_{xk}\ge 1,\forall x\in V_P,k\in V_S,\text{(24)}& {\displaystyle \sum }_{i\in V_C}{M}_i=y,\text{(25)}& M_i\in \left\{\mathrm{0,1}\right\},\hspace{1em}\forall i\in V_C,\text{(26)}& P_{xk}\in \left\{\mathrm{0,1}\right\},\hspace{1em}\forall x\in V_P,\hspace{1em}k\in V_S,\text{(27)}& H_{ij}\in \left\{\mathrm{0,1}\right\},\hspace{1em}\forall i\in V_C,\hspace{1em}j\in \left\{\mathrm{1,2},\dots ,m\right\},\end{array}$$

Formula (20) is the objective function, which consists of transportation cost and node construction cost for the intermodal logistics system. Constraint (21) reflects the capacity restrictions of underground freight handling at UD; constraint (22) ensures that the capacity of metro lines are not violated; constraint (23) stipulates that the ULS terminal on a single metro line is unique, but the UDC can be assigned to different ULS terminals based on reality; constraint (24) guarantees that UD is connected to at least one ground terminal of which ownership is unique; constraint (25) specifies that the number of metro stations that are selected as underground depots equals y; and constraints (26)–(28) are the range of relevant decision variables.

5. Solution Method

CM-ULSLAP is an NP-hard problem in which the customer layer (demand point), middle layer (GT), underground logistics layers (UD and UH), and supply layer (UDC) are integrated into a complex intermodal network. This problem is computationally intractable due to the following reasons. First, there is not a single layer of node information that can be fully invoked. Different intermodal strategies result in the conflict between customer demands and economic considerations. The evaluation procedure is hard to carry out. On the contrary, the possible paths for each OD pairs are so huge that it is computationally expensive and unsustainable to enumerate all of them for each configuration involving transport mode choice and corresponding routing transfer. Once the supply relationship and the generation strategy of candidate GT are determined, CM-ULSLAP reduces to the bi-level capacitated intermodal network design problem (BCIND). For calculation efficiency and stability, the access routes from UDCs to ULS terminals were determined in advance to cut off the supply layer. Then the entropy weighted TOPSIS method (E-TOPSIS) was applied to evaluate the demand characteristics in model (P1). Next, an exact algorithm for set coverage model is served for obtaining feasible solutions in model (P2), thereby cutting off the customer layer by replacing abundant demand points with candidate GT sets. Immune Clone Selection Algorithm (ICSA) was applied to optimize the location of UD and allocation of deterministic GT which following the minimum objective cost in model (P3). An adaptive mechanism of underground route navigation was also embedded in ICSA where fleet flows always follow the shortest path that is currently available until the certain section is blocked. Finally, the solution space of combined models has been reduced from five dimensions $M\ast N\ast K\ast \left|V_{\mathrm{C}}\right|\ast \left|V_{\mathrm{P}}\right|$ to three dimensions $M\ast \left|V_{\mathrm{C}}\right|\ast \left|V_{\mathrm{P}}^{″}\right|$ which contributes a lot to the computational efficiency. The relationships and decomposition of models are shown in Figure 2.

[figure omitted; refer to PDF]

The ICSA was inspired by the principle of clonal selection of antibodies in the biological immune system [41]. It was used to explain the basic characteristics of adaptive immune response to antigen stimulation. The ICSA has been widely used in solving problems such as combinatorial optimization, intelligent optimization, and production scheduling due to its powerful data search capabilities [42–44]. However, the convergence rate of original ICSA is slower, immune probability and cloning probability are relatively fixed, and the degree of change is relatively low when solving complex problems. In this paper, we introduced a self-adaptive mutation operation based on normal distribution to achieve a uniform and dynamic global high-level ambiguity for each antibody that meets the mutation rate for the characteristics of multiple UD feasible solutions and frequent information feedback at each stage. Probability variation enhances the randomness and stability of the search and can effectively prevent the solution from falling into local optimum. The main operation of the heuristic algorithm is shown in Figure 3.

[figure omitted; refer to PDF]

The solution steps of GT location are listed in Algorithm 1.

Algorithm 1: Set coverage model for ground terminal.

(1)

Initialization $X_j=0,Y_{ij}=0$, $V_R^{\ast }$

(2)

clear $A\left(j\right)$, $B\left(i\right)$

(3)

while $X_j=0$ do

(4)

 for $j^{\prime }\in V_R^{\ast }$, find $\left|A\left(j^{\prime }\right)\right|=\max \left\{\left|A\left(j\right)\right|\right\}$

(5)

  if true then

(6)

   let $X_{j^{\prime }}=1$, remove $j^{\prime }$ from $V_R^{\ast }$

(7)

  end if

(8)

 end for

(9)

 for $j\in A\left(j\right)$ do

(10)

  assign to $j^{\prime }$ subsequently according to length $B\left(i\right)$

(11)

  while $Q_{\max }^{j\prime }=0$ or $A\left(j^{\prime }\right)=\varnothing$ do

(12)

  end while

(13)

 for $i\in A\left(j^{\prime }\right)$ and $Y_{ij}\le 1$ do

(14)

  if $\overline{Q_i}\left(1-Y_{ij}\right)\le Q_{\max }^{j\prime }$ then

(15)

   let $Y_{i{j}^{\prime }}=1-Y_{ij}$, $Q_{\max }^{j\prime }=Q_{\max }^{j\prime }-\overline{Q_i}\left(1-Y_{ij}\right)$, $Y_{ij}=1$

(16)

   remove $i$ from $V_R^{\ast }$

(17)

  end if

(18)

  if $\overline{Q_i}\left(1-Y_{ij}\right)>Q_{\max }^{j\prime }$ then

(19)

   let $Y_{i{j}^{\prime }}=Q_{\max }^{j\prime }/\overline{Q_i}$, $Y_{ij}=Y_{ij}+Y_{i{j}^{\prime }}$, $Q_{\max }^{j\prime }=0$

(20)

  end if

(21)

check $V_R^{\ast }\ne \varnothing$

(22)

  if true then

(23)

   return to line 4

(24)

  end if

(25)

end while

Step 1: Generate Initial Population. The initial antibodies were obtained from the memory cells that were generated randomly from the feasible solution space. The coding of antibodies includes the UD open strategies and the location information of newly built GT which were fed back by E-TOPSIS. The initial antibody population was recorded as $P_0\left(\text{RC}\right)$, which was composed by RC number randomly generated antibodies.

Step 2: Diversity Evaluation. Select k individuals with the highest fitness value and t individuals with the lowest affinity from the parent population $P_n\left(\text{RC}\right)$ to compose the solution vector, where $k=\text{RC}\times 20\%$ and $t=\text{RC}\times 10\%$. Affinity was used to indicate the matching degree between antibody and antigen, namely, the optimal satisfaction between the candidate GT sets and open UD in objective (19). Lehmer mean [45] was recommended to the affinity determination. The affinity of the uth solution vector can be described as$$\begin{array}{}\text{(28)}& I_u=\frac{1}{1+L\left(u,v\right)};\forall L\left(u,v\right)=\frac{{{\displaystyle \sum }}_{u\ne v}{‖J_b\left(u\right)-J_b\left(v\right)‖}^2}{{{\displaystyle \sum }}_{u\ne v}‖J_b\left(u\right)-J_b\left(v\right)‖}.\end{array}$$

Step 3: Clone Operation. The clone ratio was calculated by the evaluation of affinity and similarity between antigen and other antibodies. The cloned population $Z_n\left(N_{\mathrm{C}}\right)$ is produced by replicating $k+t\pm d$ number of selected antibodies. Specifically, the total number of clones generated for all selected antibodies is$$\begin{array}{}\text{(29)}& N_{\text{RC}}={\displaystyle \sum }_{u=1}^{k+t\pm d}\text{round}\left(\text{RC}\left(1-I_u\right)\right).\end{array}$$

Step 4: Gene Mutation. Select a group of clones from $Z_n\left(N_{\text{RC}}\right)$ and perform Gaussian mutation on them to obtain population $S_n$. The mutation rate adopts an adaptive strategy [46] which is related to the adaptation value $f_k$ of the antibody. This mutation operation can be described as $m_j=\text{normrnd}\hspace{1em}\left(m_j,\sigma ,\mathrm{1,1}\right)$, where $m_j$ is the jth attribute of the clone, and normrnd is a normal-distribution random number with a mean of $m_j$ and a standard deviation of $\sigma$. The $\sigma$ domain of the antibody is adaptively adjusted to $\sigma =\omega \cdot I_u/f_k$ following its adaptation value and affinity.

Step 5: Immune Selection. The memory population M is consisted of a batch of antibodies with the highest $f_k$ from the set $S_n\cup Z_n$. Then select 30% k number of antibodies with the highest fitness value from M to update the equal number of antibodies with the lowest fitness value in the population $P_n$ to form the progeny population $P_{n+1}$.

Step 6: Node Search Reset. After reaching the maximum number of iterations N, check $I_y^{\left(n+1\right)}\left(i\right)$, if i exists so that $I_y^{\left(n+1\right)}\left(i\right)>\underset{t\ne y}{\max }{I}_t^{\left(n+1\right)}$. Then update the objective cost $J_b\left(i\right)$. Let $y=y+1$, returning to step 2; if not exists, step forward.

Step 7. The algorithm terminates when reaching iteration number N. Extract the most-affinity antibodies with their attributions as the preferred configuration of M-ULS.

6. Scenario Analysis

6.1. Background and Parameter Setting

The Nanjing city is located at latitude $\mathrm{N}{32}^{\circ }02^{\prime }$ and longitude $\mathrm{E}{118}^{\circ }46^{\prime }$ in East China. This city is situated on the Yangtze river bank. Based on the 2016 census, Nanjing had a population of 8.27 million. Express delivery volume of Nanjing had reached 630 million parcels in 2017 (ASKCI consulting). The total amount of social logistics is nearly 3 trillion RMB, with a 7.96% annual growth rate (Nanjing Municipal Bureau of Commerce). Three major transportation hubs, Nanjing Railway Station, Nanjing South Railway Station, and Lukou Airport, run through the city of which freight throughput capacity breaks 800,000 tons per year. Nowadays, there are 9 metro lines and 164 stations in Nanjing. The total length of metro lines is 347 kilometers, the 7th longest in the world.

For the sake of studying the performance of metro freight in the real-world case, Nanjing Metro was taken as a case to schedule the optimal configuration layout under minimum objective cost. A scenario of 2 open UDC, Line 1 and Line 2 with an activated UH and 31 stations, was involved based on the real-world location. As shown in Figure 4, there are 274 demand points distributed in a range of 201 square kilometers from Xinjiekou to Nanjing South Railway Station. Relevant data and parameter values were obtained by the Monte Carlo simulation method, as shown in Table 1. The underground freight volume of each UDC was limited to 10,000 tons per day to be in consistent with the transport capacity of underground tunnels. Cargos from two metro lines can be transferred at Xinjiekou Station where the function of UD is invalid.

[figure omitted; refer to PDF]

Table 1

Parameters of metro freight.

Parameter	Value
Freight price coefficient, $\beta$	3
Travel cost on underground network	$\text{Uniform}\left(\mathrm{1,2}\right)$
Cost of retrofitting underground depot	$\text{Uniform}\left(\mathrm{4000,6000}\right)$
Capacity of underground depot	$\text{Uniform}\left(\mathrm{3000,4500}\right)$
Capacity of ground terminal	$\text{Uniform}\left(\mathrm{2000,3000}\right)$
Freight demand of customers	$\left|\text{Normal}\left(\mathrm{50,30}\right)\right|$
Travel time on ground delivery path	$\text{Uniform}\left(\mathrm{20,200}\right)$
Maximum service radius of ground terminal	$\text{Uniform}\left(\mathrm{2,3}\right)$

The proposed solution algorithm was coded in MATLAB R2016, and all experiments were run on a desktop computer with Windows 10, Intel Core I7-7700K 4.2 GHz processor, and 32 GB of RAM. The initial population size $p=50$ and maximum iteration number $G_{\max }=100$.

6.2. Simulation Results

The E-TOPSIS method was applied to illuminate the range of system services at first. The number of qualified demand points is reduced to 228, in which the relative proximity of the most urgent area $\max {\rho }_j$ equals 0.1916, and the critical proximity $\overline{{\rho }_j}=0.1812$. Metro Line 1 saturates faster than the Line 2 with the cargo throughput of 9,747 tons per day and 7,758 tons per day, respectively. Algorithm 1 gave a set of feasible GT location demonstrated in Figure 5. The number of ground terminal that satisfies the constraints is determined as fifteen. Then the optimization results from ICSA indicate that the minimum cost of the M-ULS-integrated logistics system is 185,900 RMB per day, while the transport cost of traditional ground logistics is 198,700 RMB per day. UDs were selected to locate on the following six stations: Xinmofanmalu, Andemen, Hanzhongmen, Daxinggong, Minggugong, and Jiqingmendajie. Figure 6 depicts the optimal open strategy of UD and corresponding allocation. Minggugong Station captured the busiest UD, with the daily freight volume of 3914.9 tons, whereas Andemen Station had the lowest freight volume of 1483.9 tons with single affiliated GT. The specific node information is shown in Table 2.

[figure omitted; refer to PDF]

[figure omitted; refer to PDF]

Table 2

Summary of optimal node information.

Open UD	Freight volume (t)	Number of GT allocation	Service radius (km)	Source freight (t) (UDC 1, UDC 2)	Travel distance above (km)	Travel distance underground (km)
Xinmofanmalu	3583.5	11	2.99	(642.3, 346.6)	10.83	16.56
		13	2.86	(730.6, 328.1)	9.86	16.56
		14	2.84	(664.5, 358.6)	12.17	16.56
		15	2.66	(360.2, 154.2)	14.93	16.56
Andemen	1483.9	2	2.18	(970.1, 513.8)	9.19	28.64
Hanzhongmen	2774.3	6	1.70	(756.2, 568.8)	9.76	19.79
		10	1.92	(717.8, 731.5)	12.66	19.79
Daxinggong	2925.9	7	1.41	(855.8, 628.7)	12.07	16.54
		9	1.83	(730.2, 711.2)	9.63	16.54
Minggugong	3914.9	3	2.95	(856.5, 513.2)	15.54	16.56
		5	2.91	(255.3, 315.5)	11.28	16.56
		8	1.75	(694.4, 765.3)	10.27	16.56
		12	2.41	(196.5, 318.2)	11.48	16.56
Jiqingmendajie	2822.5	1	2.03	(622, 793.9)	11.53	26.94
		4	2.01	(696.3, 710.3)	8.83	26.94

Table 3 depicts the comparison of traffic volume from 2 UDCs to 15 GTs under traditional point-to-point distribution mode and M-ULS, where traffic volume (TV) is defined as freight OD quantity multiplied by the transport distance [47]. It can be seen that the ground traffic volume (GTV) of GT-4 has been reduced by 69.32%, which benefits the most from underground logistics. Most GTVs got a different degree of alleviation between 30% and 65%. Note that GT-15 in Figure 5 is relatively close to the freight sources, and the GTV has increased by 36.71% after adopting M-ULS. Therefore, the traditional point-to-point model is more economical for nodes near such UDCs.

Table 3

Comparison of traffic volume between point-to-point mode and M-ULS.

Open UD	GT allocation	GTV in point-to-point distribution (ton·km)	GTV in M-ULS (ton·km)	UTV in M-ULS (ton·km)	GTV alleviation rate	Underground logistics ratio (%)
Xinmofanmalu	11	5159.7	2853.7	3453.8	−44.69%	59.48
	13	4220.2	2847.4	3513.7	−32.53%	59.06
	14	5632.0	3292.4	3572.7	−41.54%	70.03
	15	1496.3	2045.7	1685.2	+36.71%	75.31
Andemen	2	10530.6	3683.3	10624.3	−65.02%	62.68
Hanzhongmen	6	8134.4	3345.0	6557.2	−58.88%	66.97
	10	9466.0	4578.8	7172.3	−51.63%	51.59
Daxinggong	7	8466.9	4614.0	6260.6	−45.51%	60.99
	9	7230.1	3481.8	5976.0	−51.84%	63.20
Minggugong	3	7584.0	5525.3	6243.7	−27.15%	57.64
	5	2296.8	2115.6	2261.8	−7.89%	75.71
	8	6532.2	3705.5	5923.6	−43.27%	57.81
	12	1634.2	1404.4	1927.2	−14.06%	61.72
Jiqingmendajie	1	10860.7	3979.8	9535.0	−63.36%	60.46
	4	10097.0	3098.2	9472.5	−69.32%	52.59

When it came to the underground traffic volume (UTV) in M-ULS, all GTs exceeded 50% of the underground logistics ratio, where the maximum reaching 75.71% at GT-5. The average above and underground mileage at Xinmofanmalu accounted for 34.03% and 65.97%, respectively. Table 4 demonstrates the overall variation of traffic volume. Benefit from the higher utilization of underground distribution, less than half of the ground delivery is required. For this part of goods that need to be delivered to the downtown areas, the original long-distance highway travel is decomposed into several subtransport processes above and underground. As a result, fewer freight vehicles are needed. The mobility on urban traffic loops will be significantly improved and some congested road sections will be effectively alleviated from traffic pressure.

Table 4

Overall variation of traffic volume.

	Overall GTV alleviation rate	Overall UTV ratio	Overall GTV ratio
Point-to-point mode M-ULS	0	0	100%
	49.09%	62.47%	37.53%

6.3. Sensitivity Analysis

Figure 7 indicates the variation of optimal cost with freight price coefficient $\beta$ and number of UD-opened $n_{\mathrm{c}}$. Note that the total cost of the intermodal system dropped slightly along with the increase of $\beta$ and then reached the same level with traditional ground logistics when $\beta$ equaled 2.52. Accordingly, there is a wide range of pricing options for underground freight transport, but the ratio of the discreet unit price of M-ULS to traditional ground logistics should be less than 1 : 2.5 in order to achieve economic advantages. On the contrary, the total cost of intermodal system ascended along with the increase of open UD, while the transport cost dropped due to the advance in network connectivity. Both trends were relatively stable. Therefore, in M-ULS-integrated logistics system, optimal cost is not sensitive to the number of UD opened. It is considerable to set up additional underground depots so that the capacity and robustness system service can be further augmented.

[figures omitted; refer to PDF]

7. Conclusions

This paper developed a multistage model with a combined algorithm based on hierarchical optimization, which aimed to deal with optimal nodes location-allocation problem in the cooperative operation of M-ULS. The applicability of the model was proven by selecting the appropriate metro lines and stations for problem solving. The critical issue of developing a M-ULS network is to establish an operational framework and optimization method for multiple nodes. Three types of nodes were defined as underground depot (UD), underground hub (UH), and ground terminal (GT), in accordance with the current state and planning of the city, logistics, and metro network. The M-ULS is in limited capacity, which has a multiple network hierarchy and multiple sources of freight. Therefore, it is necessary to screen the demand points and integrate them into the model for global optimization. A set of mixed integer programming model combined with comprehensive evaluation was developed as a more rigorous and convenient analysis approach to solve CM-ULSLAP. The model was embedded with three evaluation variables, namely, UFF, RA, and SD. A hybrid algorithm composed of E-TOPSIS, exact algorithm, and heuristic algorithm was designed to obtain the solution of the model. A real-world simulation based on Nanjing Metro was carried out to verify the effectiveness of proposed models and algorithms.

It was proved that the potential of metro stations and underground tunnels can be fully exploited to achieve higher logistics benefits. According to the results, optimized M-ULS contributed to reducing the overall ground traffic volume by 49.09% and the daily comprehensive cost by 6.44% within the service scope. Served by the new logistics system, 62.47% of the total freight was completed by ULS, and the remaining 37.53% was delivered to all demand points which were delivered by road distribution. Furthermore, compared to traditional point-to-point logistics mode, the price gap of M-ULS optimal cost has mostly remained stable, regardless of the number of UDs. The realization of these advantages neither affects the normal operation of the metro system nor further increases the amount of underground space development. The other comprehensive benefits for urban transport, environment, and society are obvious. This paper not only proposes a well-designed location-allocation model for joint operation of ULS and metro but provides a possible method for the collaborative research of different urban infrastructures. Moreover, it can provide novel ideas for the sustainable urban development.

Despite the above contributions, this study has several noteworthy limitations. First, the proposed model was exclusively applied to cross-shaped or #-shaped radial metro layout. The loop network will have some special problems that cannot be ignored, for example, the trade-off among transshipment efficiency, routing selection, and budget. Second, the proposed deterministic system in which the UHs were always opened was actually not economical. Although the new problem always starts from certainty, future research needs to explore the uncertainty problems. Robust optimization approach can be used to study the UH scheduling strategy under conditions of stochastic disruption and demand uncertainty. Additionally, research on network planning or underground space can be further extended by considering the operational conflicts between M-ULS and metro system.

Conflicts of Interest

The authors declare that they have no conflicts of interest.

Acknowledgments

This study was supported by the National Natural Science Foundation of China (71631007 and 71601095).

References

[1] Z. Chen, J. Dong, R. Ren, "Urban underground logistics system in China: opportunities or challenges?," Underground Space, vol. 2 no. 3, pp. 195-208, DOI: 10.1016/j.undsp.2017.08.002, 2017.

[2] D. Stein, R. Stein, D. Beckmann, B. Schoesser, "CargoCap: feasibility study of transporting containers through underground pipelines," Proceedings of the 4th International Symposium on Underground Freight Transportation by Capsule Pipelines and Other Tube/Tunnel Systems, .

[3] M. C. van der Heijden, A. van Harten, M. J. R. Ebben, Y. A. Saanen, E. C. Valentin, A. Verbreck, "Using simulation to design an automated underground system for transporting freight around Schiphol Airport," Interfaces, vol. 32 no. 4,DOI: 10.1287/inte.32.4.1.49, 2002.

[4] F. Cotana, F. Rossi, A. Marri, "Pipe§net: innovation in the transport through high rate small volume payloads," Proceedings of the 5th International Symposium on Underground Freight Transportation by Capsule Pipelines and Other Tube/Tunnel Systems, .

[5] S. Kashima, R. Nakamura, M. Matano, T. Taguchi, T. Shigenaga, "Study of an underground physical distribution system in a high-density, built-up area," Tunnelling and Underground Space Technology, vol. 8 no. 1, pp. 53-59, DOI: 10.1016/0886-7798(93)90137-k, 1993.

[6] H. Liu, "Feasibility of using pneumatic capsule pipelines in New York City for underground freight transport," Proceedings of the ASCE Pipeline Division Specialty Congress–Pipeline Engineering and Construction, pp. 1085-1096, .

[7] J. Visser, B. Wiegmans, B. Konings, "Review of underground logistic systems in the Netherlands: an ex-post evaluation of barriers, enablers and spin-off," Proceedings of the 5th International Symposium on Underground Freight Transportation by Capsule Pipelines and Other Tube/Tunnel Systems, .

[8] Q. Qian, "The evaluation of status quo and outlook for the underground logistics system (ULS) development in China," Proceedings of the 5th International Symposium on Underground Freight Transportation by Capsule Pipelines and Other Tube/Tunnel Systems, .

[9] K. Lewandowski, "Underground freight transport in cities, utopia or future?," Systems: Journal of Transdisciplinary Systems Science, vol. 12 no. 2, 2007.

[10] X. Zhao, B. G. Hwang, G. S. Yu, "Identifying the critical risks in underground rail international construction joint ventures: case study of Singapore," International Journal of Project Management, vol. 31 no. 4, pp. 554-566, DOI: 10.1016/j.ijproman.2012.10.014, 2013.

[11] B. G. Hwang, X. Zhao, G. S. Yu, "Risk identification and allocation in underground rail construction joint ventures: contractors’ perspective," Journal of Civil Engineering and Management, vol. 22 no. 6, pp. 758-767, DOI: 10.3846/13923730.2014.914095, 2016.

[12] M. Shan, B. G. Hwang, K. S. N. Wong, "A preliminary investigation of underground residential buildings: advantages, disadvantages, and critical risks," Tunnelling and Underground Space Technology, vol. 70, pp. 19-29, DOI: 10.1016/j.tust.2017.07.004, 2017.

[13] G. Arends, E. Boer, "Tunneling of infrastructure: from non-considered to ill considered–lessons from the Netherlands," Tunnelling and Underground Space Technology, vol. 16 no. 3, pp. 225-234, DOI: 10.1016/s0886-7798(01)00043-8, 2001.

[14] A. Dampier, M. A. Marinov, "Study of the feasibility and potential implementation of metro-based freight transportation in Newcastle upon Tyne," Urban Rail Transit, vol. 1 no. 3, pp. 164-182, DOI: 10.1007/s40864-015-0024-7, 2015.

[15] J. C. Rijsenbrij, B. A. Pielage, "The potential of metro systems for city logistics," Proceedings of the 3rd International Symposium on Underground Freight Transportation by Capsule Pipelines and Other Tube/Tunnel Systems, .

[16] J. Kikuta, T. Ito, I. Tomiyama, S. Yamamoto, T. Yamada, "New subway-integrated city logistics system," Procedia–Social and Behavioral Sciences, vol. 39, pp. 476-489, DOI: 10.1016/j.sbspro.2012.03.123, 2012.

[17] J. Kelly, M. Marinov, "Innovative Interior designs for urban freight distribution using light rail systems," Urban Rail Transit, vol. 3 no. 4, pp. 238-254, DOI: 10.1007/s40864-017-0073-1, 2017.

[18] S. J. Hong, "A proposal of green logistic system design through the subway network–focus on Seoul Metro Network," Journal of Digital Design, vol. 14 no. 3, pp. 327-336, 2014.

[19] K. Cochrane, M. J. Roorda, A. Shalaby, "Moving freight on public transit: best practices, challenges, and opportunities," International Journal of Sustainable Transportation, vol. 11 no. 2, pp. 120-132, DOI: 10.1080/15568318.2016.1197349, 2017.

[20] A. Motraghi, M. Marinov, "Analysis of urban freight by rail using event-based simulation," Simulation Modelling Practice and Theory, vol. 25, pp. 73-89, DOI: 10.1016/j.simpat.2012.02.009, 2012.

[21] E. Fatnassi, J. Chaouachi, W. Klibi, "Planning and operating a shared goods and passengers on demand rapid transit system for sustainable city-logistics," Transportation Research Part B: Methodological, vol. 81, pp. 440-460, DOI: 10.1016/j.trb.2015.07.016, 2015.

[22] M. Peng, J. Dong, R. Ren, Z. Chen, Y. Chen, "Analysis on the cooperative operation between urban underground logistics system and Metro," Journal of PLA University of Science and Technology,DOI: 10.12018/j.issn.1009-3443.20170606006/2017.09.26, 2017.

[23] B. Hörl, H. Dörr, M. Wanjek, A. Romstorfer, "METRO.FREIGHT.2020-strategies for strengthening rail infrastructure for freight transport in urban regions," Transportation Research Procedia, vol. 14, pp. 2776-2784, DOI: 10.1016/j.trpro.2016.05.478, 2016.

[24] R. Masson, A. Trentini, F. Lehuédé, N. Malhéne, O. Péton, H. Tlahig, "Optimization of a city logistics transportation system with mixed passengers and goods," EURO Journal on Transportation and Logistics, vol. 6 no. 1, pp. 81-109, DOI: 10.1007/s13676-015-0085-5, 2017.

[25] F. Zhou, G. Zhou, Z. Zhang, "The point-routing problem of urban distribution based on subway," Control and Decision, vol. 33 no. 7, 2017.

[26] D. Ambrosino, M. G. Scutellà, "Distribution network design: new problems and related models," European Journal of Operational Research, vol. 165 no. 3, pp. 610-624, DOI: 10.1016/j.ejor.2003.04.009, 2005.

[27] J. Castro, S. Nasini, F. Saldanha-Da-Gama, "A cutting-plane approach for large-scale capacitated multi-period facility location using a specialized interior-point method," Mathematical Programming, vol. 163 no. 1-2, pp. 411-444, DOI: 10.1007/s10107-016-1067-6, 2017.

[28] F. Fotuhi, N. Huynh, "A reliable multi-period intermodal freight network expansion problem," Computers & Industrial Engineering, vol. 115, pp. 138-150, DOI: 10.1016/j.cie.2017.11.007, 2018.

[29] S. A. Mirzapour, K. Y. Wong, K. Govindan, "A capacitated location-allocation model for flood disaster service operations with border crossing passages and probabilistic demand locations," Mathematical Problems in Engineering, vol. 2013,DOI: 10.1155/2013/507953, 2013.

[30] X. Li, J. Kurths, C. Gao, J. Zhang, Z. Wang, Z. Zhang, "A hybrid algorithm for estimating origin-destination flows," IEEE Access, vol. 6, pp. 677-687, 2018.

[31] R. Chauhan, T. Singh, A. Tiwari, A. Patnaik, N. S. Thakur, "Hybrid entropy—TOPSIS approach for energy performance prioritization in a rectangular channel employing impinging air jets," Energy, vol. 134, pp. 360-368, DOI: 10.1016/j.energy.2017.06.021, 2017.

[32] L. Y. Sun, C. L. Miao, L. Yang, "Ecological-economic efficiency evaluation of green technology innovation in strategic emerging industries based on entropy weighted TOPSIS method," Ecological Indicators, vol. 73, pp. 554-558, DOI: 10.1016/j.ecolind.2016.10.018, 2017.

[33] E. Nathanail, G. Adamos, M. Gogas, "A novel approach for assessing sustainable city logistics," Transportation Research Procedia, vol. 25, pp. 1036-1045, DOI: 10.1016/j.trpro.2017.05.477, 2017.

[34] J. R. Montoya-Torres, A. Muñoz-Villamizar, C. A. Vega-Mejía, "On the impact of collaborative strategies for goods delivery in city logistics," Production Planning & Control: The Management of Operations, vol. 27 no. 6, pp. 443-455, DOI: 10.1080/09537287.2016.1147092, 2016.

[35] A. V. Binsbergen, P. Bovy, "Underground urban goods distribution networks," Innovation: The European Journal of Social Sciences, vol. 13 no. 1, pp. 111-128, 2000.

[36] H. Gordijn, "Underground freight transport in the Netherlands in the next century," Journal of Economic and Social Geography, vol. 90 no. 2, pp. 234-241, DOI: 10.1111/1467-9663.00064, 1999.

[37] V. Marchau, W. Walker, R. Van Duin, "An adaptive approach to implementing innovative urban transport solutions," Transport Policy, vol. 15 no. 6, pp. 405-412, DOI: 10.1016/j.tranpol.2008.12.002, 2008.

[38] O. N. Egbunike, A. T. Potter, "Are freight pipelines a pipe dream? A critical review of the UK and European perspective," Journal of Transport Geography, vol. 19 no. 4, pp. 499-508, DOI: 10.1016/j.jtrangeo.2010.05.004, 2011.

[39] W. Winkelmans, "Redressing the balance between demand for mobility and supply of transport by means of new modes of transportation," Proceedings of the 5th International Symposium on Underground Freight Transportation by Capsule Pipelines and Other Tube/Tunnel Systems, .

[40] J. Visser, A. Binsbergen, "Underground logistical systems in cities: a visualization of the future," Proceedings of the 2nd International Symposium on Underground Freight Transportation by Capsule Pipelines and Other Tube/Tunnel Systems, .

[41] Q. Lin, Y. Ma, J. Chen, "An adaptive immune-inspired multi-objective algorithm with multiple differential evolution strategies," Information Sciences, vol. 430-431, pp. 46-64, DOI: 10.1016/j.ins.2017.11.030, 2018.

[42] D. Lizondo, S. Rodriguez, A. Will, V. Jimenez, J. Gotay, "An artificial immune network for distributed demand-side management in smart grids," Information Sciences, vol. 438, pp. 32-45, DOI: 10.1016/j.ins.2018.01.039, 2018.

[43] Z. Zhang, X. Wang, J. Lu, "Multi-objective immune genetic algorithm solving nonlinear interval-valued programming," Engineering Applications of Artificial Intelligence, vol. 67, pp. 235-245, DOI: 10.1016/j.engappai.2017.10.004, 2018.

[44] R. Shang, B. Du, H. Ma, L. Jiao, Y. Xue, R. Stolkin, "Immune clonal algorithm based on directed evolution for multi-objective capacitated arc routing problem," Applied Soft Computing, vol. 49, pp. 748-758, DOI: 10.1016/j.asoc.2016.09.005, 2016.

[45] D. H. Lehmer, "On the compounding of certain means," Journal of Mathematical Analysis and Applications, vol. 36 no. 1, pp. 183-200, DOI: 10.1016/0022-247x(71)90029-1, 1971.

[46] T. Gong, T. Fan, L. Pei, Z. Cai, "Magnetic resonance imaging-clonal selection algorithm: an intelligent adaptive enhancement of brain image with an improved immune algorithm," Engineering Applications of Artificial Intelligence, vol. 62, pp. 405-411, DOI: 10.1016/j.engappai.2016.10.004, 2017.

[47] J. Zheng, H. X. Liu, "Estimating traffic volumes for signalized intersections using connected vehicle data," Transportation Research Part C: Emerging Technologies, vol. 79, pp. 347-362, DOI: 10.1016/j.trc.2017.03.007, 2017.