Optimization of the Classification Yard Location Problem Based on Train Service Network

Abstract

This paper investigates the problems of locating yards and allocating works of train formation among yards based on the train service network. It not only involves the modifications of the scales of existing yards, either improving, or downsizing or even demolishing, but also the determinations of building new yards in a given rail network. Besides, the train service network is also taken into consideration, so that the car-hour costs that are incurred from the reclassification and accumulation operations of railcars can be optimized simultaneously. The accumulation parameter setting is symmetrical on both directions of traffic flows. A binary integer programming model is first proposed for optimizing the train service network, the solution of which is employed as a benchmark of the further integrated optimization. Based on this, a nonlinear joint optimization model is developed aiming at striking a balance between the capital investment and car-hour consumptions, with the constraints of the reclassification capacities of yards and the number of sorting tracks, and multiple logical relations among decision variables. Corresponding linearization techniques are introduced for transforming the nonlinear models into the linear ones. Finally, an exact solving approach is presented with computational results that are based on a real-world based case to illustrate the efficacy of the linear model.

1. Introduction

Classification yards (yards, terminals, marshalling stations) constitute a vital component of the rail transportation network for freight, performing the function of accumulating and reclassifying railcars. In today’s competitive environment, the design of freight service network has been concentrated on by railroads in countries around the world. As one key element of influencing the structure of the rail service network, the optimization of locations and capacities of yards in the rail system is the aim of this study, in order to improve the reliability of the service network and accelerate the transportation process of railcars.

There are 49 classification yards in China. Among them, the largest one is the Zhengzhou North marshalling station, which was once the largest yard in Asia. It has the length of 6.63 km and the width of 0.80 km, covering an area of 5.3 square kilometers. Additionally, the design capacity of the Zhengzhou North marshalling station is 36 thousand railcars per day. In accordance with the Zhengzhou railway terminal planning (2016–2030), the location of Zhengzhou North yard will be moved out due to the extension of the city. Apparently, a huge number of capitals will be invested for the relocation. Besides, the capacity of a classification yard is of great importance, since it has to be guaranteed that all of the original demands can be shipped within appropriate time. Nevertheless, if the reclassification capacities of certain yards are not sufficient, parts of demands will make detours on their trips, which might incur increased travel and sorting costs. Therefore, the location and the capacity of a classification yard are two peculiarities that can exert a significant impact on the reclassification strategies of railcars. Additionally, when the budget is limited, the investment plan of yards to be established or extended in a rail network should satisfy the budget constraint and ensure that all of the shipments can be reclassified at their perfect yards in the meanwhile.

The location optimization of yards in a rail network is actually a systemwide problem, which cannot be worked out if each yard is independently considered, since there exist interactions among yards. Thus, the efficiency of the train service network (TSN) is taken into consideration to evaluate the overall delivery efficiency of shipment.

It should be pointed out that the original demands will be disturbed if the structure of a yard network is changed and it can be calculated by the method of prediction. Either the opening or closing of a yard in the target network has an absolute influence on the traffic volumes of demands. The revenue in aspects of serving the demands will also change in this situation, which can be employed as a feedback of optimizing the classification yard location problem (CYLP) and be used for modifying the optimization result of CYLP. The simultaneous optimization of revenue in aspects of serving the demands and CYLP is intractable. Only the method of dynamic optimization is qualified to take this revenue into account. However, the optimization of the location problem (LP) of classification yards is to satisfy the current transportation demands. Additionally, the efficiency of the TSN is calculated based on the current original demands instead of the changed predict demands. Thereby, the change of demand and the revenue in aspects of serving the demands are not considered.

When there exist a huge number of yards in the rail network, the number of investment combinations for establishing or extending yards is exponential. For example, the number of combinations will be 1024 if there are five yards in a simple rail network and for each yard there are four investment plans; if there are ten yards and four investment plans for each yard, the number of combinations will grow to 1,048,576. Furthermore, each alternative investment combination is corresponding to a reclassification strategy of shipments. Hence, the optimization of classification yard location is such a complicated and practical problem that requires being addressed immediately.

The facility LP has been studied by plenty of scholars with their sustained enthusiasm over the years. Alfred initially introduced the location theory [1], with the aim of finding an optimal location of a single warehouse and minimizing the total travel costs between the warehouse and customers. This work has then inspired the subsequent researchers (e.g., Isard [2]; O’Kelly [3]). Surveys of relevant literatures can be found in Brandeau and Chiu [4], who have investigated over 50 representative problems in location research, and Revelle et al. [5], in which the modeling approach is divided into four categories: analytic, continuous, network, and discrete location models, and the LP is classified into median and plant LPs, and center and covering problems. It is well known that the LP has very popular applications in transportation (e.g., Klincewicz [6]; Nickel et al. [7]; Gelareh and Nickel [8]), communication (e.g., Schneider and Zastrow [9]; Pirkul, [10]), and other network design problem that requires efficient distributions. With special regards to rail transportation, Mansfield and Wein [11] presented a pioneering model for the classification yards’ LP in rail, which is applied to serve a railroad company in choosing appropriate locations for new facilities. However, this model was primarily constructed to solve the problem of yard location, without the consideration of reclassification strategy of railcars, which led to a tremendous reclassification cost. Lin et al. [12] investigated the decision of the classification yard location from a pre-defined alternative set for a new yard, which was also predetermined to be built in a rail network, the decision of the investment plan for an existing yard was also explored. A bi-level multi-project programming model was then constructed with the aim of minimizing the total cost of investment, reclassification, and accumulation based on the design of binary variables of train service, flow combination, and investment plan. Their study also considered the influence of reconstruction or building of an existing or new yard on original demand by introducing a constant coefficient. However, this coefficient cannot describe the change of the yard network, since it is a constant and irrelevant to the original demands. Besides, the demolishment of a yard was not considered in their study. Moreover, the solving approach of this nonlinear programming model was not introduced, which makes the feasibility of their proposed model unknown. Jeong et al. [13] addressed a hub-and-spoke network problem for railroad freight, which was usually favored by air freight. Their final solution revealed potential hub locations. However, this kind of network structure is not actually employed widely in rail. Lee et al. [14] combined the influence of consolidating of flows and developed an optimal marshalling yard location model that was based on the economies of scale.

In general, the LP is involved with locating facilities and allocating lower level nodes (e.g., demand nodes) to higher level nodes (e.g. hubs). Related studies can be referred to Eben-Chaime et al. [15], Alumur et al. [16]. The classification yards actually have no concept of hierarchical but only the difference of scales and automation due to the peculiarity of railway network. Consequently, the allocation problem in rail LP is referred to the arrangement of train services in the rail network and the assignment of reclassification workload (e.g., Liu et al., [17]). The most relevant studies to our problem are those jointly given the locations of yards and the allocation of reclassification workloads. Lin [18] elaborated the bi-level programming model based on the previous works. A new binary variable was separated from the original investment variable that was designed in Lin et al. [12] in order to describe the situation of employing no investment plan on the existing yard. Whereas, in this present work, the yards that remain unchanged are classified into one set, which have already been predetermined. The constraints considered in the lower-level model in Lin [18] included the reclassification capacity and track number constraints, but the constraint of the number of blocks that can be formed in a classification yard is not taken into account. Besides, similar to Lin et al. [12], the solving approach is not mentioned neither. In the work of Liu et al. [19], the planning horizon was divided into two period and the original demands were changed in the second period according to the investment plan that was obtained in period. The proposed linear integer model is not exactly a dynamic programming, since only two timestamps are actually designed in their work. However, it can provide inspiration for the dynamic optimization of CYLP based on the prediction of demands. Besides, only the improvement of an existing yard or building of a new yard is involved in Liu et al. [19], the downsizing and closure are not permitted. Nevertheless, the situations of both downsizing and demolishment of yards are taken into account in our research. Furthermore, the commercial software was employed directly and a rail network with only two candidate yards for improvement was devised for the illustration of feasibility and effectiveness of their model due to the linear peculiarity of their proposed model. Moreover, Liu et al. [17] presented a simulated annealing algorithm to tackle the multi-classification-yard location (MCYL) problem that was employed on a large-scale China railway network.

There is a rich body of researches that deal with the design of the TSN. As one of the prior works on the train service network design problem (TSNDP), Assad [20] first presented a hierarchical taxonomy of modelling issues in rail from a holistic perspective, including strategic, tactical, and operational decisions. From the angles of network flows and combinatorial optimization, typical models dealing with car routing and train makeup were presented. The resource allocation problem in rail was also discussed, which was defined as the assignment of cars to block and blocks to trains. The basic concepts of rail service network design were described in this research, which can provide an instructional and inspirational model for other scholars. Crainic et al. [21] examined the problems of routing shipments, scheduling train services, and allocating classification works between yard pairs. A nonlinear mixed integer multicommodity flow model was proposed. In this model, the limitation of the capacities of train services were taken into consideration, while that of the reclassification capacities of yards is more practical. Furthermore, based on these previous researches, Haghani [22] combined the distribution problem of empty cars with train routing and dispatching problems. With regard to the recent researches, Lin et al. [23] tackled the large-scale train connection services problem with simulated annealing algorithm with a bi-level programming model. This model can determine the optimal sequence of train services assigned to each shipment and the frequencies of services. Yaghini et al. [24] investigated the train formation planning problem, as a special case of service network design problems, corresponding to the determination of the assembled trains, including the definition of assignments, frequencies, sizes, and profiles of shipments and traction. An improved local branching approach was proposed to deal with the proposed 0–1 integer programming model.

The contributions of this paper can be summarized, as follows:

• The location-allocation problem of classification yards in a rail network is investigated combined with the TSNDP. This combination and the interaction between the CYLP and TSNDP are elaborated on the basis of several small-scale artificial cases, so that the significance and necessity of addressing CYLP can be clarified.
• Not only the scenario of improving the scale of an existing yard is considered, but also that of downsizing or demolishing a yard is taken into account. A set of the decision variables are designed in order to determine whether an original existing yard or a new yard to-be-built exists after optimization in the given network.
• Two nonlinear programming models are proposed aiming at striking a balance between economic costs and car-hour consumptions in rail. Moreover, two linearization approaches are presented to transform the nonlinear models into the linear ones.

The rest of this paper is organized, as follows: Section 2 presents a comprehensive literature review related to CYLP. Section 3 elaborates the CYLP based on a small-scale artificial network. Section 4 constructs several nonlinear programming models for CYLP, proposes corresponding linearization approaches and conducts a numerical example and Section 5 gives the conclusion of this paper.

2. The Classification Yard Location Problem

Five small-scale artificial networks are endowed in Figure 1, Figure 2, Figure 3, Figure 4 and Figure 5 to more specifically describe the CYLP. Four existing classification yards Y1, Y2, Y3, and Y4 (concentric circles) and several stations around these yards (small circles, which are used to express the structure of the rail network) are included in each network. Train services (blocks) are expressed as blue dotted and orange dashed arcs with arrows.

In the first four networks, there is a positive traffic flow between each origin-destination (OD) pair. It is worth pointing out that the OD demand that is generated inside a yard is not under consideration, because it involves the optimization within the yard instead of between the yards. Therefore, the volumes of these inside flows are regarded as zero. Usually, there are two specific shipping patterns for the commodities that are delivered among yards. Either they will be assigned to train services so that these shipments can be shipped directly to their destinations, or they will be reclassified at the intermediate yards. When it comes to the second pattern, portions of reclassification capacity of a yard will be consumed. Additionally, this situation of capacity consumption, of course, should be taken into consideration. While, in the following examples, it is assumed that the capacities of all the yards are large enough to reclassify the illustrative flows, which is, the constraint of reclassification capacity is not considered for now.

Another important thing is the limitation of the number of the sorting tracks in a yard. A relatively relaxed restriction of sorting tracks is that the total number of them should be larger than the number of occupied tracks by railcars. For example, if one sorting track of a yard can handle 200 railcars and the reclassification workload of this yard is 450 railcars per day. Subsequently, the number of the needed sorting tracks would be three tracks (the number of tracks must be integer). Comparatively, a rigorous restriction of sorting tracks is that the number of them cannot be larger than the number of train services (blocks) generated at one yard. We assume that it is acceptable if one block occupies one or several sorting tracks, while not if more than one block occupies one sorting track. In fact, the rigorous restriction is more reasonable and practical, and it is considered in the examples in Figure 1, Figure 2, Figure 3, Figure 4 and Figure 5. The numbers of available sorting tracks (excluding these for local trains) of yards are presented above them. For example, yard Y1 has three sorting tracks, Y2 five, Y3 five, and Y4 eight before any modifications.

Notice that the TSNs that are exhibited in Figure 1, Figure 2, Figure 3, Figure 4 and Figure 5 are all artificial without optimization. These service network might not be the optimal, but we just use these train services to illustrate the CYLP. In Figure 1, the numbers of occupied tracks are three, five, five, and three for yard Y1, Y2, Y3, and Y4, respectively. Moreover, there is a total of eight train services provided and the utilizations of sorting tracks reach to 100% for Y1, Y2, and Y3, and five tracks remain for Y4. Figure 1 offers a perfect example for the TSNDP, since it is the basic network before optimizing with the CYLP simultaneously.

According to the shipping strategy that is given below the network, there is a total of 12 traffic flows in the network and four of them need to be reclassified at the intermediate yards Y2 and Y3. The total reclassification times are four. Actually, this shipping strategy only presents one potential transportation pattern of commodities. There are still several alternative choices of shipping for shipments. Taking shipment $f_{14}$ as an example, it is first assigned to the train service form Y1 to Y3 and then placed on the train service from Y3 to Y4 to finish the delivery according to the given strategy presented in Figure 1. However, it is also possible that this flow is shipped by the train services from Y1 to Y2, Y2 to Y3, and then Y3 to Y4. Which shipping strategy to choose depends on the incurred accumulation delay (cost) and the reclassification cost. For example, if the traffic volume of the shipment $f_{14}$ is 150.00 cars per day, and the accumulation parameter of yard Y1 and Y3 is 12.0 hours (notice that the accumulation parameters might be influenced by the directions of the traffic flows. For example, if a specific bulk commodity is generally shipped from north to south, then this bulk car flow with the direction from north to south is easier to be accumulated than the opposite direction. Nevertheless, such influence is not considered in this work and we assume that the accumulation parameter is symmetrical on both directions of flows). The reclassification cost per railcar is 4.2 and 4.3 hours for Y2 and Y3. Subsequently, the accumulation delay (the calculation of accumulation delay can be referred to Lin et al. [23]) of this flow will be 12.0 $\times$ 50 $\times$ 2 = 1200.0 car-hour and the reclassification delay will be (4.2 + 4.3) $\times$ 150.00 = 1275 car-hour. Since the accumulation delay of $f_{14}$ is larger than the reclassification delay, shipment $f_{14}$ should assigned to the train services from Y1 to Y3 and Y3 to Y4 instead of being reclassified twice if it is assigned to the train services from Y1 to Y2, Y2 to Y3, and Y3 to Y4.

Figure 2 describes the scenario of improving yard Y2 by adding an extra sorting track. In fact, if a yard is planned to be improved or downsized, not only the number of sorting tracks will be increased, but the reclassification capacity will also be enhanced. Nevertheless, in the illustrative examples, the reclassification capacities are ignored for the sake of simplification. After the improvement, yard Y2 has six sorting tracks and one more train service is provided from Y2 to Y4. In this scenario, yard Y4 has four occupied sorting tracks out of eight. The shipment $f_{24}$ can then be shipped by the train service from Y2 to Y4 directly instead of being reclassified at Y3. The total reclassification times reduce to three in this situation.

Similar to the scenario in Figure 3, yard Y1 is improved by increasing two more train services from Y1 to Y4 and Y4 to Y1, respectively. In this case, yard Y4 has five occupied tracks out of eight. The shipments $f_{14}$ and $f_{41}$ are assigned to the new train services, respectively. The total reclassification times are two in this case. Scenario in Figure 4 combines the two improvement plans of Y1 and Y2 and the total reclassification time reduces to 1 in this case. It can obviously illustrate that, the more the train service provided, the less the reclassification times are. However, providing more train services might not be the most economical and car-hour-saving approach to deliver the shipments due to the accumulation delay that is incurred by forming these train services. It has to be simultaneously optimized with the reclassification operations.

Figure 5 explores the scenario of building a new yard Y5 on the basis of the original network and there is no improvement plans for other existing yards. The need for building a new yard will occur if the existing yards cannot undertake the current reclassification tasks. Generally, newly-built yards are of high-efficiency and they can share more workloads than the existing yards. It can also offer alternative shipping choices for the OD demands. Besides, if a new yard is established in a network, it will attract more original demands from the surrounding stations. For example, the construction of the new yard Y5 will generate two more OD demands, $f_{15}$ and $f_{54}$. Because the constraints of the sorting track numbers of the existing yards, only two new train services are formed in the network, from Y1 to Y5 and Y5 to Y4.

In the five scenarios that are depicted in Figure 1 to Figure 5, only the improvements of the existing classification yards are considered, while the downsizing or even demolishment are not, which could be likely in the real world. Therefore, in the modelling process, the possibilities of downsizing and closure of existing yards are included. Moreover, the restrictions mentioned above, such as the reclassification capacities and the number of sorting tracks, can be regarded as the practical constraints. There are also other critical logical constraints, for example, if the demand $f_{14}$ in the Figure 1 prefers to be reclassified at yard Y1, then the train service from yard Y1 to Y3 must be provided. Such logical constraints also need to be taken into consideration. The detailed mathematical model is introduced later in Section 3.

The reclassification problem is actually very complicated due to the difficulties of routing shipments in a rail network. For example, in a line rail network, if five yards are included in a line rail network, then there will be 150 combinations. Furthermore, if there are eight marshalling yards, there will be 1.3 billion combinations for routing shipments (see Lin et al., [23]). The number of combinations of railcars is an explosive growth with the number of yards in the network, not to mention the integrated optimization with the investment plans of yards. Therefore, the CYLP is worth studying and it has practical values.

3. Mathematical Model for the CYLP

This section introduces two nonlinear 0–1 integer programming models, including several continuous auxiliary variables for the TSNDP and CYLP, respectively. Corresponding linearization approaches are proposed to transform the nonlinear models into the linear ones, so that an exact solving approach can be employed.

3.1. Notations

In this study, the reconstruction plans for the existing yards and the building plan for the new yards in a rail network are both considered. Moreover, in each reconstruction plan, not only the alternative plan in which the reclassification capacity is extended and the sorting tracks’ number is increased for an existing yard is taken into account, but also the potential plan in which the sorting capacity is reduced and the sorting tracks’ number is decreased for an existing yard is considered. Additionally, even each existing yard might be demolished in accordance with their function and performance in a rail network. Of course, the corresponding capital investments of different reconstruction strategies are different. Besides, the scales of the new yards to-be-built have multiple choices. Therefore, all of the yards in a rail network can basically be divided into three classes: yards that stay the current states and are without any changes; existing yards that need to be either improved, or reduced, or demolished; or, potential new yards that might be built, based on which the subsequent optimization and discussion are proposed and presented.

Table 1. Notations employed in the models.

Notation	Description
Sets
$V^{\mathrm{All}}$	The set of all the yards, including yards to-be-built and reconstructed;
$V^{\mathrm{Stay}}$	The set of the existing yards that stay unchanged, that is, without any reconstructions, in a rail network;
$V^{\mathrm{Modified}}$	The set of the existing yards that might require reconstructions or demolitions in a rail network;
$V^{\mathrm{New}}$	The set of the yards which are possible to be built in a rail network;
$plan(k)$	The investment plan of yard $k$;
$\vartheta (i,j)$	The set of yards through which a flow from $i$ to $j$ may pass excluding yard $i$ and yard $j$.
Parameters
$c_i$	The accumulation parameter of yard $i$; (hour)
$m_{ij}$	The size of a train dispatched from yard $i$ to $j$; (cars)
${\tau }_k^{}$	The classification cost for each railcar at yard $k$; (hours per railcar)
${\tau }_k^{p-\mathrm{New}}$	The classification cost for each railcar at the new yard $k$ in accordance with the investment plan $p$ if yard $k$ is determined to be built after optimization; (hours per railcar)
${\tau }_k^{p-\mathrm{Modified}}$	The classification cost for each railcar at the reconstructed yard $k$ in accordance with the investment plan $p$; (hours per railcar)
$C_k^{Total}$	The original total reclassification capacity of the existing yard $k$; (cars per day)
$\beta$	The proportion of reserved reclassification capacity for the local trains at the existing yard $k$;
${\rho }_k^{p-\mathrm{New}}$	The reclassification capacity of the new yard $k$ according to the investment plan $p$ if yard $p$ is optimized to be built; (cars per day)
${\rho }_k^{p-\mathrm{Modified}}$	The changed (increased or reduced) reclassification capacity (excluding for the local trains) of yard $k$ according to the investment plan $p$ if yard $k$ needs to be reconstructed after optimizing; (cars per day)
$T_k^{\mathrm{Total}}$	The original number of sorting tracks of the existing yard $k$; (tracks)
$\theta$	The proportion of reserved sorting tracks for the local railcars at the existing yard $k$;
${\sigma }_k^{p-\mathrm{New}}$	The number of sorting tracks of the new yard to be built in line with the investment plan $p$ after optimization; (tracks)
${\sigma }_k^{p-\mathrm{Modified}}$	The changed (increased or reduced) number of sorting tracks of yard $k$ in line with the investment plan $p$ after optimization; (tracks)
$M$	An extremely huge number, which is used to express the logical relations among variables;
$I_k^p$	The cost of investment plan $p$ of yard $k$; (CNY)
$\alpha$	The converting coefficient that unifies conflicting units of costs; (car-hour per CNY)
$\omega$	The budget of investment; (CNY)
$N_{ij}$	The original demand which originates at yard $i$ and is destined for yard $j$.
Decision variables
$y_k^p$	Investment variable; it takes value one if plan $p$ is selected for yard k; otherwise, it is zero;
${\varsigma }_{ij}$	Train service variable; it takes the value of one if the train service from yard $i$ to j is provided; otherwise, it is zero;
$r_{ij}^k$	Car flow variable; it takes value one if railcars destinating for yard j areconsolidated into train service i→k at yard $i$; Otherwise, it is zero.
$x_i$	Yard variable; it takes value one if yard $i$ is existing or newly established; otherwise, it is zero.
Auxiliary variables
fij	The actual traffic flow from yard $i$ to j including the original demand from i to j and the reclassified flows sorted at yard i from other yards (notice that the origins of the reclassified flows could not be i);
Dij	The service flow from yard $i$ to j, which is used to indicate the traffic volume shipped by train service from i to j.

first lists the notations employed in our models.

3.2. The Nonlinear 0-1 Programming Model for CYLP

Before introducing the mathematical model for rail CYLP, a general model (TSNDP) was employed to assign (allocate) the reclassification workload to rail yards or tracks is presented. The aim of this model is to provide a benchmark before optimizing.

(nonlinear TSNDP, M1)

$$minimize{\displaystyle \sum _{i\ne j\hspace{0.33em},i,j\in V^{\mathrm{All}}-V^{\mathrm{New}}}{c}_i{m}_{ij}{\varsigma }_{ij}}+{\displaystyle \sum _{k\ne i,\hspace{0.33em}k\ne j,\hspace{0.33em}k\in V^{\mathrm{All}}-V^{\mathrm{New}}}{\displaystyle \sum _{i\ne j\hspace{0.33em},i,j\in V^{All}-V^{\mathrm{New}}}{f}_{ij}{r}_{ij}^k{\tau }_k}}$$

(1)

$$subject to {\varsigma }_{ij}+{\displaystyle \sum _{k\ne i,\hspace{0.33em}k\ne j,\hspace{0.33em}k\in \vartheta (i,j)}{r}_{ij}^k}=1 \forall i\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}}-V^{\mathrm{New}}$$

(2)

$$r_{ij}^k\le {\varsigma }_{ik}\hspace{1em}\forall k\ne i,\hspace{0.33em}k\ne j,\hspace{0.33em}i\ne j,\hspace{0.33em}i,j,k\in V^{\mathrm{All}}-V^{\mathrm{New}}$$

(3)

$${\displaystyle \sum _{i\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}}-V^{\mathrm{New}}}{f}_{ij}}{r}_{ij}^k\le (1-\beta )C_k^{Total} \forall k\ne i,\hspace{0.33em}k\ne j,\hspace{0.33em}k\in V^{\mathrm{All}}-V^{\mathrm{New}}$$

(4)

$${\displaystyle \sum _{j\in V^{\mathrm{All}}}\psi (D_{ij})}\le (1-\theta )T_i^{Total} \forall i\in V^{\mathrm{All}}-V^{\mathrm{New}}$$

(5)

$$\sum _{j\in V^{\mathrm{All}}}{\varsigma }_{ij}(1-\le )T_i^{Total}\hspace{1em}\forall i\in V^{\mathrm{All}}-V^{\mathrm{New}}$$

(6)

$$D_{ij}=f_{ij}{\varsigma }_{ij}+{\displaystyle \sum _{t\ne i,t\ne j,t\in V^{\mathrm{All}}-V^{\mathrm{New}}}{f}_{it}}{r}_{it}^j \forall i\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}}-V^{\mathrm{New}}$$

(7)

$$f_{ij}=N_{ij}+{\displaystyle \sum _{s\ne i,s\ne j,s\in V^{\mathrm{All}}-V^{\mathrm{New}}}{f}_{sj}}{r}_{sj}^i \forall i\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}}-V^{\mathrm{New}}$$

(8)

$${\varsigma }_{ij}\in \left\{0,1\right\}\hspace{1em}\forall i\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}}-V^{\mathrm{New}}$$

(9)

$$r_{ij}^k\in \left\{0,1\right\}\hspace{1em}\forall i\ne j,\hspace{0.33em}k\ne i,\hspace{0.33em}k\ne j\hspace{0.33em}i,j,k\in V^{\mathrm{All}}-V^{\mathrm{New}}$$

(10)

$$f^{ij},D^{ij}\ge 0\hspace{1em}\forall i\ne j,i,j\in V^{\mathrm{All}}-V^{\mathrm{New}}$$

(11)

The objective function is to minimize the total cost that is incurred by the reclassification operations and accumulation processes of railcars. The basic information of shipping strategies of traffic flows, the reclassification workloads of yards and sorting tracks occupied by railcars, and the construction of TSN is obtainable by solving this model. The solution of M1 actually can be regarded as a benchmark that is intended to be compared with the solution achieved under the case of modifying the configuration of yards or building new yards in a rail network. The CYLP model is then able to be constructed in the basis of TSNDP (M1). The specific meanings of the constraints in M1 are introduced in the subsequently constructed CYLP model since parts of the constraints of TSNDP and CYLP model have the analogous implications.

The motivation of the CYLP model is to minimize the total capital investment of yards and operation costs from the perspective of railroads. Based on the model M1, the nonlinear 0–1 programming mathematical model for the CYLP is presented, as follows:

(nonlinear CYLP, M2)

$$\begin{array}{ll}minimize& {\displaystyle \sum _{k\in V^{\mathrm{Modified}}\cup V^{\mathrm{New}}}{\displaystyle \sum _{p\in plan(k)}\alpha I_k^p{y}_k^p}}+{\displaystyle \sum _{i\ne j\hspace{0.33em},i,j\in V^{All}}{c}_i{m}_{ij}{\varsigma }_{ij}}{x}_i{x}_j+{\displaystyle \sum _{k\ne i,\hspace{0.33em}k\ne j,\hspace{0.33em}k\in V^{\mathrm{Stay}}}{\displaystyle \sum _{i\ne j\hspace{0.33em},i,j\in V^{All}}{f}_{ij}{r}_{ij}^k{\tau }_k}}+\\ & {\displaystyle \sum _{k\ne i,\hspace{0.33em}k\ne j,\hspace{0.33em}k\in V^{\mathrm{New}}}{\displaystyle \sum _{i\ne j\hspace{0.33em},i,j\in V^{\mathrm{All}}}{\displaystyle \sum _{p\in plan(k)}{f}_{ij}{r}_{ij}^k{y}_k^p{\tau }_k^{p-\mathrm{New}}}}}+{\displaystyle \sum _{k\ne i,\hspace{0.33em}k\ne j,\hspace{0.33em}k\in V^{\mathrm{Modified}}}{\displaystyle \sum _{i\ne j\hspace{0.33em},i,j\in V^{\mathrm{All}}}{\displaystyle \sum _{p\in plan(k)}{f}_{ij}{r}_{ij}^k{y}_k^p{\tau }_k^{p-\mathrm{Modified}}}}}\end{array}$$

(12)

$$subject\hspace{0.33em}to\hspace{1em}{x}_i=1\forall i\in V^{\mathrm{Stay}}$$

(13)

$${\varsigma }_{ij}\le x_i\hspace{1em}\forall i\ne j,\hspace{0.33em}i\in V^{\mathrm{New}}\cup V^{\mathrm{Modified}},\hspace{0.33em}j\in V^{\mathrm{All}}$$

(14)

$${\varsigma }_{ij}\le x_j\hspace{1em}\forall i\ne j,\hspace{0.33em}j\in V^{\mathrm{New}}\cup V^{\mathrm{Modified}},\hspace{0.33em}i\in V^{\mathrm{All}}$$

(15)

$${\varsigma }_{ij}+{\displaystyle \sum _{k\ne i,\hspace{0.33em}k\ne j,\hspace{0.33em}k\in \vartheta (i,j)}{r}_{ij}^k}=1\hspace{1em}\forall i\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}}$$

(16)

$$r_{ij}^k\le {\varsigma }_{ik}\hspace{1em}\forall k\ne i,\hspace{0.33em}k\ne j,\hspace{0.33em}i\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}},k\in V^{\mathrm{All}}$$

(17)

$${\displaystyle \sum _{i\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}}}{f}_{ij}{r}_{ij}^k}\le (1-\beta )C_k^{Total}\hspace{1em}\forall k\ne i,\hspace{0.33em}k\ne j,\hspace{0.33em}k\in V^{\mathrm{Stay}}$$

(18)

$${\displaystyle \sum _{i\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}}}{f}_{ij}{r}_{ij}^k}\le {\displaystyle \sum _{p\in plan(k)}{\rho }_k^{p-\mathrm{New}}{y}_k^p}\hspace{1em}\forall k\ne i,\hspace{0.33em}k\ne j,\hspace{0.33em}k\in V^{\mathrm{New}}$$

(19)

$${\displaystyle \sum _{p\in plan(k)}{y}_k^p}\le x_k\hspace{1em}\forall k\in V^{\mathrm{New}}$$

(20)

$${\displaystyle \sum _{i\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}}}{f}_{ij}{r}_{ij}^k}\le (1-\beta )C_k^{Total}+{\displaystyle \sum _{p\in plan(k)}{\rho }_k^{p-\mathrm{Modified}}{y}_k^p}\hspace{1em}\forall k\ne i,\hspace{0.33em}k\ne j,\hspace{0.33em}k\in V^{\mathrm{Modified}}$$

(21)

$${\displaystyle \sum _{p\in plan(k)}{y}_k^p}\le x_k\hspace{1em}\forall k\in V^{\mathrm{Modified}}$$

(22)

$${\displaystyle \sum _{j\in V^{\mathrm{All}}}\psi (D_{ij})}\le (1-\theta )T_i^{Total}\hspace{1em}\forall i\in V^{\mathrm{Stay}}$$

(23)

$${\displaystyle \sum _{j\in V^{\mathrm{All}}}{\varsigma }_{ij}}\le (1-\theta )T_i^{Total}\hspace{1em}\forall i\in V^{\mathrm{Stay}}$$

(24)

$${\displaystyle \sum _{j\in V^{\mathrm{All}}}\psi (D_{ij})}\le {\displaystyle \sum _{p\in plan(i)}{y}_i^p}{\sigma }_i^{p-\mathrm{New}}\hspace{1em}\forall i\in V^{\mathrm{New}}$$

(25)

$${\displaystyle \sum _{j\in V^{\mathrm{All}}}{\varsigma }_{ij}}\le {\displaystyle \sum _{p\in plan(i)}{y}_i^p}{\sigma }_i^{p-\mathrm{New}}\hspace{1em}\forall i\in V^{\mathrm{New}}$$

(26)

$${\displaystyle \sum _{j\in V^{\mathrm{All}}}\psi (D_{ij})}\le (1-\theta )T_i^{Total}+{\displaystyle \sum _{p\in plan(i)}{y}_i^p{\sigma }_i^{p-\mathrm{Modified}}}\hspace{1em}\forall i\in V^{\mathrm{Modified}}$$

(27)

$${\displaystyle \sum _{j\in V^{\mathrm{All}}}{\varsigma }_{ij}}\le (1-\theta )T_i^{Total}+{\displaystyle \sum _{p\in plan(i)}{y}_i^p{\sigma }_i^{p-\mathrm{Modified}}}\hspace{1em}\forall i\in V^{\mathrm{Modified}}$$

(28)

$${\displaystyle \sum _{k\in V^{\mathrm{Modified}}\cup V^{\mathrm{New}}}\hspace{1em}{\displaystyle \sum _{p\in plan(k)}{I}_k^p{y}_k^p\le }}\hspace{0.33em}\omega$$

(29)

$$D_{ij}=f_{ij}{\varsigma }_{ij}+{\displaystyle \sum _{t\ne i,t\ne j,t\in V^{\mathrm{All}}}{f}_{it}}{r}_{it}^j\hspace{1em}\forall i\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}}$$

(30)

$$f_{ij}=N_{ij}+{\displaystyle \sum _{s\ne i,s\ne j,s\in V^{\mathrm{All}}}{f}_{sj}}{r}_{sj}^i\hspace{1em}\forall i\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}}$$

(31)

$$y_k^p\in \left\{0,1\right\}\hspace{1em}\forall k\in V^{\mathrm{New}}\cup V^{\mathrm{Modified}},\hspace{0.33em}p\in plan(k)$$

(32)

$${\varsigma }_{ij}\in \left\{0,1\right\}\hspace{1em}\forall i\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}}$$

(33)

$$r_{ij}^k\in \left\{0,1\right\}\hspace{1em}\forall i\ne j,\hspace{0.33em}k\ne i,\hspace{0.33em}k\ne j\hspace{0.33em}i,j,k\in V^{\mathrm{All}}$$

(34)

$$x_i\in \left\{0,1\right\}\hspace{1em}\forall i\in V^{\mathrm{All}}$$

(35)

$$f_{ij},D_{ij}\ge 0\hspace{1em}\forall i\ne j,i,j\in V^{\mathrm{All}}$$

(36)

The objective function (12) aims at minimizing the total cost of capital investments of yards, accumulation, and reclassification delays of railcars. Since the optimization of rail yards’ location is in the position of railroad, the units of capital investments of yards are converted into car-hour from CNY. Besides, the reclassification delays are calculated in three terms, according to different properties of yards.

Equation (13) is the existing yard constraint. It ensures that the existing yards that have no requirements to be reconstructed or demolished are bound to be remaining.

Equations (14) and (15) are the train service constraints. They guarantee that the train services can only be provided if both yard $i$ and $j$ are included in the rail network after optimizing.

Equation (16) is the exclusive condition. It restricts that only one single transportation pattern can be selected for each shipment from yard $i$ to $j$, which is, this traffic flow is either to be shipped directly from yard $i$ to $j$ without any reclassification operations, or to be transported to the front yard $k$ first on its itinerary, and then consolidated with other flows at yard $k$.

Equation (17) ensures that only if the train service from yard $i$ to $k$ is provided, shipment originating from yard $i$ and destinated at yard $j$ can be reclassified at yard $k$ on its itinerary.

Equations (18), (19), and (21) are the reclassification constraints of yards. They guarantee that the total volume of the reclassified traffic flow cannot exceed the capacity of each yard. Notice that a portion of the reclassification capacity should be reserved for the local trains. Usually, it accounts for 15 percent of total reclassification capacity in China. Moreover, Equations (20) and (22) guarantee that, at most, one of the investment options is selected for a yard $k$, provided that it becomes part of the network.

Equations (23)–(28) are the sorting track constraints of yards. They ensure that both the number of the sorting tracks occupied by all train services and that of the formed blocks are less than that of the available sorting tracks of yards. Similarly, part of the sorting tracks is supposed to be reserved for the local railcars. Generally, it consumes 15 percent of the whole sorting tracks, the same as the proportion of reserved reclassification capacity. The number of sorting tracks occupied by service flows at yard $i$ can be calculated by:

$${\displaystyle \sum _{j\in V^{\mathrm{All}}}\psi (D_{ij})}={\displaystyle \sum _{j\in V^{\mathrm{All}}}\left(\left[D_{ij}/\omega \right]+1\right)}$$

(37)

where $\omega$ denotes the maximum number of railcars that can be stored on a sorting track. The parameter $\omega$ usually takes the value of 200 cars, according to the design standard. Besides, it is unpractical and uneconomical to place more than one destination-oriented train service on a reclassification track (referred to Assad [20]). One outbound train can pick up one or multiple train services from different reclassification tracks according to its “take-list”, but cannot pick up multiple train services from one sorting track at a given yard. Consequently, the number of computational tracks should round up to an integer.

Equation (29) describes the limitation of capital investments. The total investment costs should be larger than the deterministic budget.

Equation (30) indicates the compositions of service flows. For a service flow $D_{ij}$ from yard $i$ to $j$, it includes not only the actual traffic $f_{ij}{\varsigma }_{ij}$, which is delivered by bypass trains directly from yard $i$ to $j$ without intermediate reclassifications, but also other actual flows $f_{it}^j$ originating from yard $i$, but first reclassified at yard $j$ on their itineraries. Equation (31) demonstrates that the total flow of traffic from yard $i$ to $j$ is the number of cars either originating before $i$ or destined at $j$. Figure 6 can illustrate these two types of traffic flows.

In Figure 6, there are eight yards (concentric rings) and seven traffic flows (black arcs with arrows) in the illustrative rail network. Among these flows, $N_{45}$ is the original demand from yard Y4 to Y5, $f_{15}$, $f_{25}$, and $f_{35}$ are the traffic flows reclassified at Y4, $f_{46}$, $f_{47}$, and $f_{48}$ are the traffic flow originating from Y4 and classified at Y5. The pre-given train services (only train services among Y1, Y2, Y3, Y4, and Y5 are known, expressed as the green dotted lines with arrows, while train services among Y5, Y6, Y7, and Y8 are unclear in this case) are from yard Y1 to Y4, Y2 to Y4, Y3 to Y4, and Y4 to Y5, say, $y_{14},\hspace{0.33em}{y}_{24},\hspace{0.33em}{y}_{34}\hspace{0.33em}and\hspace{0.33em}{y}_{45}=1$. Since the reclassification cost at yard Y4 is calculated in accordance with the total flows handled at this yard, all of the traffic flows that are placed on the outbound trains at Y4 should be counted. Consequently, the actual flow $f_{45}$ should be the sum of traffic flows $N_{45}$, $f_{15}$, $f_{25}$, and $f_{35}$, that is, $f_{45}=N_{45}+f_{15}+f_{25}+f_{35}$. The service flow from Y4 to Y5 is made up of all the traffics that are assigned to train service from Y4 to Y5. Therefore, the service flow $D_{45}$ includes traffic flows $f_{45}$, $f_{46}$, $f_{47}$, and $f_{48}$, which is, $D_{45}=f_{45}+f_{46}+f_{47}+f_{48}$.

Equations (32)–(36) are the variable domain constraints.

The above model can be employed to determine: (a) the location of yards and the allocation of reclassification workload; (b) the construction of service network in a given rail network; and, (c) the shipping strategies of traffic flows. This model can straightforwardly demonstrate the real-world restrictions and costs at classification yards and simultaneously optimize the location of yards. Nevertheless, model M2 is a nonlinear programming that cannot be directly solved through high-efficiency linear programming solvers, such as Gurobi. Therefore, several processing techniques are employed in the next subsection to linearize the model.

3.3. The Linearization Techniques

The objective function (1), and the constraints (4), (5), (7), (8) are nonlinear, as presented in M1. In model M2, it can be observed that the objective function (12), constraints (18), (19), (21), (23), (25), (27), (30), and (31) are all the nonlinear ones. Basically, there are two forms of nonlinear terms in these equations: the product of two binary variables, i.e., $x_i\times x_j,\hspace{0.33em}{x}_i,\hspace{0.33em}{x}_j\in \left\{0,1\right\}$ and the product of binary and continuous variables, say, $x\times f,\hspace{0.33em}x\in \left\{0,1\right\},\hspace{0.33em}f\ge 0$. Consequently, two corresponding linearization approaches are discussed. Here, we assume that there are multiple binary variables $x_1,\hspace{0.33em}{x}_2,\cdots ,\hspace{0.33em}{x}_n$, and let the $X=x_1\times x_2\times \cdots \times x_n$. Subsequently, the expression of the binary variable $X$ can be transformed to the following form (a similar linearization technique can also be referred to Liu et al. [19], Zhao and Lin, [25]):

$${\displaystyle \sum _{i=1}^n{x}_i}-(n-1)\le X\le \frac{1}{n}{\displaystyle \sum _{i=1}^n{x}_i}$$

(38)

There are three scenarios that are involved in Equation (38): (a) all of the values of binary variables $x_i\hspace{0.33em}(i=1,2,\cdots ,n)$ are equal to one; (b) all of the binary variables $x_i\hspace{0.33em}(i=1,2,\cdots ,n)$ take the value of zero; and, (c) parts of the values of binary variables $x_i\hspace{0.33em}(i=1,2,\cdots ,n)$ equal to zero while others one. For scenario (a), the value of the binary variable $X$ is equal to one, which is in accordance with Equation (38): $1=n-(n-1)\le X\le n/n=1$. As for scenario (b), the binary variable $X$ should equal to zero, which can be concluded from Equation (38) as well ($-(n-1)\le 0\le X\le 0/n=0$). With regard to scenario (c), the binary variable $X$ still equals zero, even if portions of $x_i\hspace{0.33em}(i=1,2,\cdots ,n)$ are positive. Let $a$ ($a<n$) denote the total number of variables $x_i\hspace{0.33em}(i=1,2,\cdots ,n)$ with the value of one. Subsequently, Equation (38) can be specified as: $a-(n-1)\le 0\le X\le a/n<1$. Consequently, Equation (38) involves all of the possible values of the product of binary variables.

Subsequently, another linearization technique can be demonstrated by introducing a binary variable $y$ and a continuous variable $f$ and let $Y$ denote the product of $y$ and $f$. The formulation of the continuous variable $Y$ can be re-expressed by:

$$f-M(1-y)\le Y\le f+M(1-y)$$

(39)

$$-My\le Y\le My$$

(40)

When the binary variable $y$ equals zero, the continuous variable $Y$ will also be zero. Under this circumstance, Equations (39) and (40) will be re-written as $f-M\le 0\le Y\le M\le f+M$ and $0\le Y\le 0$, which induces that $Y$ equals to zero, say $y=0,\hspace{0.33em}Y=0$. As the binary variable, $y$ takes the value of one, the value of the continuous variable $Y$ will be equal to that of the variable $f$. In this case, Equations (39) and (40) can be specified to: $f\le Y\le f$ and $-M\le Y\le M$, which is, $y=1,\hspace{0.33em}Y=f$.

On the basis of the above linearization techniques, several binary and continuous variables are designed to transform the nonlinear expressions in M1 and M2:

$$f_{ij}^k=f_{ij}{r}_{ij}^k\hspace{1em}\forall k\ne i,\hspace{0.33em}k\ne j,\hspace{0.33em}i\ne j,\hspace{0.33em}i,j,k\in V^{\mathrm{All}}-V^{\mathrm{New}}$$

(41)

$$F_{ij}=f_{ij}{\varsigma }_{ij}\hspace{1em}\forall i\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}}-V^{\mathrm{New}}$$

(42)

$${\xi }_{ij}={\varsigma }_{ij}{x}_i{x}_j\hspace{1em}\forall i\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}}$$

(43)

$$f_{ij}^k=f_{ij}{r}_{ij}^k\hspace{1em}\forall k\ne i,\hspace{0.33em}k\ne j,\hspace{0.33em}i\ne j,\hspace{0.33em}i,j,k\in V^{\mathrm{All}}$$

(44)

$$F_{ij}=f_{ij}{\varsigma }_{ij}\hspace{1em}\forall i\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}}$$

(45)

$${\mathcal{F}}_{ij}^k=f_{ij}^k{x}_k\hspace{1em}\forall k\ne i,\hspace{0.33em}k\ne j,\hspace{0.33em}i\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}},\hspace{0.33em}k\in V^{\mathrm{All}}$$

(46)

$${\mathbb{F}}_{ij}^{kp}={\mathcal{F}}_{ij}^k{y}_k^P\hspace{1em}\forall k\ne i,\hspace{0.33em}k\ne j,\hspace{0.33em}i\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}},\hspace{0.33em}k\in V^{\mathrm{New}}\cup V^{\mathrm{Modified}}$$

(47)

The nonlinear 0–1 programming model M1 can be reformulated as:

(linear TSNDP, M3)

$$minimize{\displaystyle \sum _{i\ne j\hspace{0.33em},i,j\in V^{\mathrm{All}}-V^{New}}{c}_i{m}_{ij}{\varsigma }_{ij}}+{\displaystyle \sum _{k\ne i,\hspace{0.33em}k\ne j,\hspace{0.33em}k\in V^{\mathrm{All}}-V^{New}}{\displaystyle \sum _{i\ne j\hspace{0.33em},i,j\in V^{\mathrm{All}}-V^{\mathrm{New}}}{f}_{ij}^k{\tau }_k}}$$

(48)

$$subject\hspace{0.33em}to\hspace{1em}\left(2\right), \left(3\right), \left(5\right), \left(6\right), \left(9\right) ~ \left(11\right)$$

$$f_{ij}-M(1-r_{ij}^k)\le f_{ij}^k\le f_{ij}+M(1-r_{ij}^k)\hspace{1em}\forall k\ne i,\hspace{0.33em}k\ne j,\hspace{0.33em}i\ne j,\hspace{0.33em}i,j,k\in V^{\mathrm{All}}-V^{\mathrm{New}}$$

(49)

$$-M{r}_{ij}^k\le f_{ij}^k\le M{r}_{ij}^k\hspace{1em}\forall k\ne i,\hspace{0.33em}k\ne j,\hspace{0.33em}i\ne j,\hspace{0.33em}i,j,k\in V^{\mathrm{All}}-V^{\mathrm{New}}$$

(50)

$$f_{ij}-M(1-{\varsigma }_{ij})\le F_{ij}\le f_{ij}+M(1-{\varsigma }_{ij})\hspace{1em}\forall i\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}}-V^{\mathrm{New}}$$

(51)

$$-M{\varsigma }_{ij}\le F_{ij}\le M{\varsigma }_{ij}\hspace{1em}\forall i\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}}-V^{\mathrm{New}}$$

(52)

$${\displaystyle \sum _{i\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}}}{f}_{ij}^k}\le (1-\beta )C_k^{Total}\hspace{1em}\forall k\ne i,\hspace{0.33em}k\ne j,\hspace{0.33em}k\in V^{\mathrm{All}}$$

(53)

$$D_{ij}=F_{ij}+{\displaystyle \sum _{t\ne i,t\ne j,t\in V^{\mathrm{All}}-V^{\mathrm{New}}}{f}_{it}^j}\hspace{1em}\forall i\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}}-V^{\mathrm{New}}$$

(54)

$$f_{ij}=N_{ij}+{\displaystyle \sum _{s\ne i,s\ne j,s\in V^{\mathrm{All}}-V^{\mathrm{New}}}{f}_{sj}^i}\hspace{1em}\forall i\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}}-V^{\mathrm{New}}$$

(55)

$$f_{ij}^k\ge 0\hspace{1em}\forall k\ne i,\hspace{0.33em}k\ne j,\hspace{0.33em}i\ne j,i,j\in V^{\mathrm{All}}-V^{\mathrm{New}}$$

(56)

$$F_{ij}\ge 0\hspace{1em}\forall i\ne j,i,j\in V^{\mathrm{All}}-V^{\mathrm{New}}$$

(57)

The nonlinear 0–1 programming model M2 can be re-formulated as:

(linear CYLP, M4)

$$\begin{array}{l}minimize\hspace{0.33em}{\displaystyle \sum _{k\in V^{\mathrm{Modified}}\cup V^{\mathrm{New}}}{\displaystyle \sum _{p\in plan(k)}\alpha I_k^p{y}_k^p}}+{\displaystyle \sum _{i\ne j\hspace{0.33em},i,j\in V^{\mathrm{All}}}{c}_i{m}_{ij}{\xi }_{ij}}+{\displaystyle \sum _{k\ne i,\hspace{0.33em}k\ne j,\hspace{0.33em}k\in V^{\mathrm{Stay}}}{\displaystyle \sum _{i\ne j\hspace{0.33em},i,j\in V^{\mathrm{All}}}{f}_{ij}^k}{\tau }_k}+\\ \hspace{0.33em}\hspace{0.33em} {\displaystyle \sum _{k\ne i,\hspace{0.33em}k\ne j,\hspace{0.33em}k\in V^{\mathrm{New}}}{\displaystyle \sum _{i\ne j\hspace{0.33em},i,j\in V^{\mathrm{All}}}{\displaystyle \sum _{p\in plan}{\mathbb{F}}_{ij}^{kp}{\tau }_k^{p-\mathrm{New}}}}}+{\displaystyle \sum _{k\ne i,\hspace{0.33em}k\ne j,\hspace{0.33em}k\in V^{\mathrm{Modified}}}{\displaystyle \sum _{i\ne j\hspace{0.33em},i,j\in V^{\mathrm{All}}}{\displaystyle \sum _{p\in plan}{\mathbb{F}}_{ij}^{kp}{\tau }_k^{p-\mathrm{Modified}}}}}\end{array}$$

(58)

$$subject\hspace{0.33em}to\hspace{1em}\left(13\right) ~ \left(17\right), \left(20\right), \left(22\right) ~ \left(30\right), \left(33\right) ~ \left(38\right)$$

$${\varsigma }_{ij}+x_i+x_j-2\le {\xi }_{ij}\le {\varsigma }_{ij}\hspace{1em}\forall i\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}}$$

(59)

$${\varsigma }_{ij}+x_i+x_j-2\le {\xi }_{ij}\le x_i\hspace{1em}\forall i\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}}$$

(60)

$${\varsigma }_{ij}+x_i+x_j-2\le {\xi }_{ij}\le x_j\hspace{1em}\forall i\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}}$$

(61)

$$f_{ij}-M(1-r_{ij}^k)\le f_{ij}^k\le f_{ij}+M(1-r_{ij}^k)\hspace{1em}\forall i\ne j,\hspace{0.33em}k\ne i,\hspace{0.33em}k\ne j,\hspace{0.33em}i,j,k\in V^{\mathrm{All}}$$

(62)

$$-M{r}_{ij}^k\le f_{ij}^k\le M{r}_{ij}^k\hspace{1em}\forall i\ne j,\hspace{0.33em}k\ne i,\hspace{0.33em}k\ne j,\hspace{0.33em}i,j,k\in V^{\mathrm{All}}$$

(63)

$$f_{ij}-M(1-{\varsigma }_{ij})\le F_{ij}\le f_{ij}+M(1-{\varsigma }_{ij})\hspace{1em}\forall i\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}}$$

(64)

$$-M{\varsigma }_{ij}\le F_{ij}\le M{\varsigma }_{ij}\hspace{1em}\forall i\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}}$$

(65)

$$f_{ij}^k-M(1-x_k)\le {\mathcal{F}}_{ij}^k\le f_{ij}^k+M(1-x_k)\hspace{1em}\forall i\ne j,\hspace{0.33em}k\ne i,\hspace{0.33em}k\ne j,\hspace{0.33em}i,j,k\in V^{\mathrm{All}}$$

(66)

$$-M{x}_k\le {\mathcal{F}}_{ij}^k\le M{x}_k\hspace{1em}\forall i\ne j,\hspace{0.33em}k\ne i,\hspace{0.33em}k\ne j,\hspace{0.33em}i,j,k\in V^{\mathrm{All}}$$

(67)

$${\mathcal{F}}_{ij}^k-M(1-y_k^P)\le {\mathbb{F}}_{ij}^{kp}\le {\mathcal{F}}_{ij}^k+M(1-y_k^P)\hspace{1em}\forall i\ne j,\hspace{0.33em}k\ne i,\hspace{0.33em}k\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}},\hspace{0.33em}k\in V^{\mathrm{New}}\cup V^{\mathrm{Modified}},\hspace{0.33em}p\in plan(k)$$

(68)

$$\begin{array}{c}-M{y}_k^P\le {\mathbb{F}}_{ij}^{kp}\le M{y}_k^P\\ \forall i\ne j,\hspace{0.33em}k\ne i,\hspace{0.33em}k\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}},\hspace{0.33em}k\in V^{\mathrm{New}}\cup V^{\mathrm{Modified}},\hspace{0.33em}p\in plan(k)\end{array}$$

(69)

$${\displaystyle \sum _{i\ne j\hspace{0.33em},i,j\in V^{\mathrm{All}}}{f}_{ij}^k}\le (1-\beta )C_k^{Total}\hspace{1em}\forall k\ne i,\hspace{0.33em}k\ne j,\hspace{0.33em}k\in V^{\mathrm{Modified}}$$

(70)

$${\displaystyle \sum _{i\ne j\hspace{0.33em},i,j\in V^{\mathrm{All}}}{f}_{ij}^k}\le {\displaystyle \sum _{p\in plan(k)}{\rho }_k^{p-\mathrm{New}}{y}_k^p}\hspace{1em}\forall k\in V^{\mathrm{New}}$$

(71)

$${\displaystyle \sum _{i\ne j\hspace{0.33em},i,j\in V^{\mathrm{All}}}{f}_{ij}^k}\le (1-\beta )C_k^{Total}+{\displaystyle \sum _{p\in plan(k)}{\rho }_k^{p-Modified}{y}_k^p}\hspace{1em}\forall k\in V^{\mathrm{Modified}}$$

(72)

$$D_{ij}=F_{ij}+{\displaystyle \sum _{t\in V^{\mathrm{All}}}{f}_{it}^j}\hspace{1em}\forall i\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}}$$

(73)

$$f_{ij}=N_{ij}+{\displaystyle \sum _{s\in V^{\mathrm{All}}}{f}_{sj}^i}\hspace{1em}\forall i\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}}$$

(74)

$${\xi }_{ij}\in \left\{0,1\right\}\hspace{1em}\forall i\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}}$$

(75)

$$F_{ij}\ge 0\hspace{1em}\forall i\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}},\hspace{0.33em}k\in V^{\mathrm{All}}$$

(76)

$${\mathcal{F}}_{ij}^k,f_{ij}^k\ge 0\hspace{1em}\forall i\ne j,\hspace{0.33em}i,j\in V^{\mathrm{All}},\hspace{0.33em}k\in V^{\mathrm{All}}$$

(77)

$${\mathbb{F}}_{ij}^{kp}\ge 0\hspace{1em}\forall i\ne j,\hspace{0.33em}i,j\in V^{All},\hspace{0.33em}k\in V^{\mathrm{New}}\cup V^{\mathrm{Modified}}$$

(78)

The mathematical model M2 is a strict linear programming model, which can be directly solved by the standard commercial solvers, i.e. CPLEX or Gurobi. A numerical example is presented in the next section to illustrate the efficacies of the linear programming models M3 and M4.

4. Numerical Example

In this section, a real-world based example is conducted in order to illustrate the efficacies of the proposed models and solving approaches.

4.1. The Input Data and Computational Results for M3

A real-world based rail network in China is constructed to text the models that are presented above and illustrate the efficacies of them furthermore. Figure 7 depicts the structure of this network. The original network consists of five existing yards (Wuwei South, Yingshuiqiao, Baotou West, Taiyuan North, and Qingdao West, marked as Y1, Y2, Y3, Y5, and Y6 by concentric rings, respectively, in Figure 7), and one potential yard (Yanan, Y4, whether to be built or not needs optimization), and several stations (small circles around yards, which are added into the network for the completion of the network).

Table 2. The volume of railcars for each OD demand (cars per day).

Yard	Y1	Y2	Y3	Y4	Y5	Y6
Y1	0.00	79.51	84.23	—	54.72	89.52
Y2	96.03	0.00	13.06	—	84.85	46.44
Y3	18.86	44.89	0.00	—	92.26	17.28
Y4	—	—	—	0.00	—	—
Y5	62.18	118.67	95.94	—	0.00	24.94
Y6	71.44	68.73	69.87	—	87.29	0.00

lists the corresponding information of traffic flows.

Notice that yard Y4 has not been built yet in this network. Consequently, there are yet no traffic flows arriving at or departing from Y4 and it cannot be served as an intermediate classification yard, which is, any one of shipments cannot be reclassified at this yard. Therefore, only the shipments delivered between two yards are considered, while the shipments generated inside yards are not.

The original demands are randomly generated with the units of railcars per day. There is a total of 20 positive OD demands in this network. Since the shipping strategies of traffics between yards are our focus, the shipments delivered within a yard are not considered. Consequently, the original demands within yards are set to zero. The total and the average OD demands are 1,320.71 and 52.83 cars per day, respectively, as listed in Table 2. Among these positive traffics, the largest demand is 118.67 from Y5 to Y2, while the smallest one is 13.06 from Y2 to Y3. The size of trains dispatched from each yard is set to 50 cars. The parameter $M$ is set to 1,000,000.

Table 3. The parameter settings.

Yard	${\mathit{c}}_{\mathit{i}}$	${\mathit{\tau }}_{\mathit{k}}$	$(1-\mathit{\beta }){\mathit{C}}_{\mathit{k}}^{\mathit{T}\mathit{o}\mathit{t}\mathit{a}\mathit{l}}$	$(1-\mathit{\theta }){\mathit{T}}_{\mathit{i}}^{\mathit{T}\mathit{o}\mathit{t}\mathit{a}\mathit{l}}$
Y1	11.0	4.1	450	15
Y2	10.9	3.8	330	15
Y3	11.5	4.2	640	12
Y5	11.8	3.9	560	7
Y6	10.5	4.4	400	11
Y1	11.0	4.1	450	15

lists other settings of necessary parameters.

The linear programming model M3 are solved by the commercial software Gurobi 8.0.1 on a 3.10 GHz Intel (R) Core (TM) i5-7276U CPU computer with 4.0 GB of RAM. A total of 94 continuous variables and 64 binary variables are generated. The best value of the objective function is 9602.34 car-hour. Figure 8 depicts the computational results.

The green dotted arcs indicate the train services provided between yards, and Y1 ~ Y6 refers to the existed or to-be-built classification yards in the network in Figure 8. Notice that a train service is bound to be provided for each adjacent yard pair. Besides, it is interesting to find out that if a train service is provided from one yard to another, it cannot be certain that the train service must be provided in the opposite direction. The train service from Y6 to Y1 has the longest distance, as shown in Figure 8. Though the travel distance is not taken into consideration in this research, the optimization goal of simultaneously minimizing the reclassification and accumulation costs restricts that the train service from Y6 to Y1 has to be provided while the train service from Y1 to Y6 has not. Moreover, Figure 8 also depicts the shipping strategies of each shipment. Taking $f_{15}$ as an example, it will be placed on the train service from Y1 to Y3 first, and then from Y3 to Y5 heading for the destination. The shipping strategies of shipments delivered between two yards are not demonstrated, since they are certainly transported directly without any intermediate reclassifications. It can be observed that each shipment is reclassified not more than once, which illustrates that the transportation of shipments achieved by the solving of model M3 is considerably efficient.

Table 4. The shipping strategies of shipments.

No.	Origin	Destination	${\mathit{N}}_{\mathit{i}\mathit{j}}$	${\mathit{f}}_{\mathit{i}\mathit{j}}$	Train Service Sequence
1	Y1	Y2	79.51	79.51	(Y1→Y2)
2	Y1	Y3	84.23	84.23	(Y1→Y2) & (Y3→Y2)
3	Y1	Y5	54.72	54.72	(Y1→Y2) & (Y2→Y5)
4	Y1	Y6	89.52	89.52	(Y1→Y2) & (Y2→Y6)
5	Y2	Y1	96.03	186.33	(Y2→Y1)
6	Y2	Y3	13.06	97.29	(Y2→Y3)
7	Y2	Y5	84.85	139.57	(Y2→Y5)
8	Y2	Y6	46.44	135.96	(Y2→Y6)
9	Y3	Y1	18.86	18.86	(Y3→Y2) & (Y2→Y1)
10	Y3	Y2	44.89	44.89	(Y3→Y2)
11	Y3	Y5	92.26	92.26	(Y3→Y5)
12	Y3	Y6	17.28	17.28	(Y3→Y6)
13	Y5	Y1	62.18	62.18	(Y5→Y1)
14	Y5	Y2	118.67	118.67	(Y5→Y2)
15	Y5	Y3	95.94	165.81	(Y5→Y3)
16	Y5	Y6	24.94	24.94	(Y5→Y6)
17	Y6	Y1	71.44	71.44	(Y6→Y2) & (Y2→Y1)
18	Y6	Y2	68.73	68.73	(Y6→Y2)
19	Y6	Y3	69.87	69.87	(Y6→Y5) & (Y5→Y3)
20	Y6	Y5	87.29	87.29	(Y6→Y5)

lists the detailed reclassification strategies obtained from M3.

The detailed reclassification workloads of yards and sorting tracks are presented in the second and the fourth columns of

Table 5. The yards and tracks workloads.

Yard	Reclassified Railcars	Capacity Utilization	Tracks Occupied	Remaining Tracks
Y1	307.98	68.44%	9	6
Y2	318.77	96.60%	10	5
Y3	173.29	27.08%	4	8
Y5	371.6	66.36%	7	0
Y6	297.33	74.33%	8	3
Ave.	293.79	66.56%	8	5

, respectively.

The Yard Y3 has the heaviest reclassification task, which is 448 cars per day, while Y2 has the highest capacity utilization, which is 96.60%, according to Table 5. With regards to the utilization of sorting tracks, there is no more free tracks for yard Y5. Either the reclassification capacities or the sorting tracks has been used sufficiently, which will lead to a low efficiency of classification even the whole system crash when more delivery demands occur. Consequently, the scenario of extended OD demands is considered, and we provide several modification plans for these existing yards and building plans for new yards to re-allocate the workload in this rail network.

4.2. The Input Data and Computational Results for M4

In this real-world based example, we assume that yards Y2, Y3, and Y5 are required to be reconstructed (extended or reduced), and the necessary of building the new yard Y4 needs to be determined. Three plans are provided for each of these yards, whether the extension or reduction for the existing yards and three construction scales for the new yard. The accumulation parameter of Y4 takes the value of 12.0.

Table 6. The additional parameter settings for M4.

	Y2	Y3	Y5		Y4
${\tau }_k^{1-\mathrm{Modified}}$	3.4	3.8	3.5	${\tau }_k^{1-\mathrm{New}}$	4.0
${\tau }_k^{2-\mathrm{Modified}}$	3.6	4.0	3.7	${\tau }_k^{2-\mathrm{New}}$	4.5
${\tau }_k^{3-\mathrm{Modified}}$	4.0	4.4	4.1	${\tau }_k^{3-\mathrm{New}}$	3.5
${\rho }_k^{1-\mathrm{Modified}}$	230	8	80	${\rho }_k^{1-\mathrm{New}}$	400
${\rho }_k^{2-\mathrm{Modified}}$	91	−92	−20	${\rho }_k^{2-\mathrm{New}}$	300
${\rho }_k^{3-\mathrm{Modified}}$	−109	−292	−220	${\rho }_k^{3-\mathrm{New}}$	450
${\sigma }_k^{1-\mathrm{Modified}}$	2	1	6	${\sigma }_k^{1-\mathrm{New}}$	20
${\sigma }_k^{2-\mathrm{Modified}}$	1	−2	2	${\sigma }_k^{2-\mathrm{New}}$	15
${\sigma }_k^{3-\mathrm{Modified}}$	−4	−6	−2	${\sigma }_k^{3-\mathrm{New}}$	25
$I_k^1$	0.2	0.1	0.2	$I_k^1$	0.4
$I_k^2$	0.1	0.2	0.3	$I_k^2$	0.2
$I_k^3$	0.4	0.6	0.2	$I_k^3$	0.3

lists other detailed parameters.

The linear programming CYLP model M4 is also solved by Gurobi 8.0.1. There are 69 continuous variables and 64 binary variables in ·total and the minimum value of the objective function is 13,747.40 car-hour.

Gurobi 8.0.1 also solves the linear programming CYLP model M4. There are 69 continuous variables and 64 binary variables in total and the minimum value of the objective function is 13,747.40 car-hour. In accordance with the computational results, none of the existing yards is demolished and the new yard Y4 should be built, so that the workloads of yards and tracks can be almost evenly allocated. Figure 9 depicts the new shipping strategies of shipments achieved by M4.

The red dashed lines in Figure 9 demonstrate the new added train services after optimization. When compared with the train services network in Figure 8, the service network in Figure 9 is much more complicated. Another nine train services are also added in the network, which indicates that the more shipments are placed on one train service and transported directly from their origins to the corresponding destinations, instead of being reclassified at the intermediate stations on their itineraries. The variation of TSN induces the changing of shipping strategies of portions of traffic flows. For example, shipment $f_{13}$ was placed on train service from Y1 to Y3 without any intermediate reclassifications, whereas it is now shipped by train service from Y1 to Y2 and Y2 to Y3. The shipments delivered between two adjacent yards are not depicted in Figure 9, since these flows are bound to be directly delivered without reclassifications. There are three new shipments originates at or heading for newly built yard Y4, which are $f_{14}$, $f_{41}$, $f_{46}$, and $f_{64}$. Among these four flows, the shipment $f_{46}$ is placed on the train service from Y4 to Y6, while the shipment $f_{14}$ and $f_{41}$ has to be reclassified at Y2 on its itinerary and $f_{64}$ is reclassified at Y5. Detailed information of shipping strategies of all the shipments is listed in

Table 7. The shipping strategies of shipments.

No.	Ori.	Des.	${\mathit{N}}_{\mathit{i}\mathit{j}}$	${\mathit{f}}_{\mathit{i}\mathit{j}}$	Train Service Sequence
1	Y1	Y2	79.51	79.51	(Y1→Y2)
2	Y1	Y3	84.23	84.23	(Y1→Y2) & (Y2→Y3)
3	Y1	Y4	78.54	78.54	(Y1→Y2) & (Y2→Y4)
4	Y1	Y5	54.72	54.72	(Y1→Y2) & (Y2→Y5)
5	Y1	Y6	89.52	89.52	(Y1→Y6)
6	Y2	Y1	96.03	111.59	(Y2→Y1)
7	Y2	Y3	13.06	97.29	(Y2→Y3)
8	Y2	Y4	78.55	157.09	(Y2→Y4)
9	Y2	Y5	84.85	139.57	(Y2→Y5)
10	Y2	Y6	46.44	46.44	(Y2→Y5) & (Y5→Y6)
11	Y3	Y1	18.86	18.86	(Y3→Y1)
12	Y3	Y2	44.89	44.89	(Y3→Y2)
13	Y3	Y4	85.45	85.45	(Y3→Y4)
14	Y3	Y5	92.26	92.26	(Y3→Y5)
15	Y3	Y6	17.28	17.28	(Y3→Y6)
16	Y4	Y1	15.56	15.56	(Y4→Y2) & (Y2→Y1)
17	Y4	Y2	25.43	25.43	(Y4→Y2)
18	Y4	Y3	35.44	35.44	(Y4→Y3)
19	Y4	Y5	28.74	28.74	(Y4→Y5)
20	Y4	Y6	78.55	78.55	(Y4→Y6)
21	Y5	Y1	62.18	62.18	(Y5→Y1)
22	Y5	Y2	118.67	187.40	(Y5→Y2)
23	Y5	Y3	95.94	165.81	(Y5→Y3)
24	Y5	Y4	87.45	123.23	(Y5→Y4)
25	Y5	Y6	24.94	71.38	(Y5→Y6)
26	Y6	Y1	71.44	71.44	(Y6→Y1)
27	Y6	Y2	68.73	68.73	(Y6→Y5) & (Y5→Y2)
28	Y6	Y3	69.87	69.87	(Y6→Y5) & (Y5→Y3)
29	Y6	Y4	35.78	35.78	(Y6→Y5) & (Y5→Y4)
30	Y6	Y5	87.29	87.29	(Y6→Y5)

, below.

The CYLP model aims at simultaneously optimizing the capital investment of reconstruction or construction, the accumulation and reclassification cost of railcars. The conversion coefficient $\alpha$ is set to 1000, which can ensure the magnitudes of CNY and car-hour costs are at the same level, in order to convert CNY to car-hour. The investment plans for yards obtained is that yards Y2 and Y5 adopt the first investment plans ($y_2^1=1$,$y_5^1=1$) and yard Y4 adopts the second investment plan ($y_4^2=1$).

Table 8. The modified workloads of yards and tracks.

Yard	Reclassified Railcars	Capacity Utilization	Tracks Occupied	Remaining Tracks
Y1	386.52	85.89%	11	4
Y2	551.98	98.57%	10	7
Y3	258.74	40.43%	5	7
Y4	183.72	61.24%	4	11
Y5	610.00	95.31%	11	2
Y6	333.11	83.28%	9	2
Ave.	387.35	77.45%	9	6

lists the reclassification workloads and the number of used tracks of yards. Yard Y2 has the highest capacity utilization among these yards, which is 98.57% and one more sorting track remains after improving yards Y2 and Y5, and building Y4.

5. Conclusions

In this paper, the CYLP in rail is defined as a particular LP with the rail characteristics. Not only the decision of the new yards’ locations is involved in CYLP, but also the allocation of train services and workloads of yards are included. A nonlinear 0–1 programming model is proposed to tackle this problem in the aim of minimizing the sum of capital investment, accumulation cost, and reclassification cost, with the constraints of reclassification capacity, the number of sorting tracks, and the investment budget. Several linearization techniques are employed to transform the nonlinear model into the linear one. The computational results can provide multiple useful information for decision makers. For example, the optimized TSNs can be employed to decide which bypass train should be provided and the origin and corresponding destination of it. The shipping strategies of traffic flows can also be concluded from the results. Besides, the critical decision of the reconstruction plans of the existing yards and the building plans of the new yards can be determined. A real-world case is presented with five existing yards and one potential yard to-be-built to test the effectiveness and validity of our model. The result indicates that our model can be served as an aid in determining the yard location and allocation of train service and reclassification workloads.