An Energy-Efficient Timetable Optimization Approach in a Bi-DirectionUrban Rail Transit Line: A Mixed-Integer Linear Programming Model

The quick growth of energy consumption in urban rail transit has drawn much attention due to the pressure of both operational cost and environmental responsibilities. In this paper, the timetable is optimized with respect to the system cost of urban rail transit, which pays more attention to energy consumption. Firstly, we propose a Mixed-Integer Non-Linear Programming (MINLP) model including the non-linear objective and constraints. The objective and constraints are linearized for an easier process of solution. Then, a Mixed-Integer Linear Programming (MILP) model is employed, which is solved using the commercial solver Gurobi. Furthermore, from the viewpoint of system cost, we present an alternative objective to optimize the total operational cost. Real Automatic Fare Collection (AFC) data from the Changping Line of Beijing urban rail transit is applied to validate the model in the case study. The results show that the designed timetable could achieve about a 35% energy reduction compared with the maximum energy consumption and a 6.6% cost saving compared with the maximum system cost.


Introduction
Compared to other modes of urban transport, the metro system has a larger transport capacity, high security, and is much more energy efficient [1]. Therefore, it has received rapid development around the world in recent years. For example, 36 cities in China operated the urban rail transit over 4800 km in 2018. According to the Web of China Rail Transit, there will be more than 50 cities operating urban rail transit over the next few years. In 2020, the total mileage of urban rail transit in China will reach 6000 km. However, Since the fast extension of the urban rail transit network, there has been a rapid corresponding increase of energy consumption in the urban rail transit system. This has attracted wide attention from both the environmental sectors and academic research [2]. In Beijing, as shown in Figure 1, from 2010-2017, the energy consumption of the urban rail transit system had increased two times to over 1.8 billion kWh in 2017 according to the data from Beijing urban rail transit. The enormous energy consumption generates not only huge operation cost, but also tough environmental responsibilities. The electricity cost accounts for over one-third of the total operation cost of Beijing urban rail transit, which has undertaken an energy-saving task from the Beijing government for the carbon emission requirement and sustainable development. Similarly, in the European Union (EU), metropolitan transportation is responsible for about 25% of the total CO 2 emissions in the transport system [1]. Moreover, there are detailed targets for the EU to reduce Timetabling is done after the upper management stage (line planning) that determines the routes and stop schedules of the trains. The line planning is at a more macroscopic level. The timetabling is a more detailed implementation, which may be adjusted annually. A well-designed timetable will have essential influences on the scheduling staff and rolling stock, which are simultaneously determined for daily operation. The basic consideration for the timetable problem is passenger flow. The passenger will profit from a well-designed train timetable with elaborate arrival times, departure times, dwell times, and running times. However, the operation cost for the operator is of equal importance. The energy consumption is not only one main cost, but also entails responsibility for environmental protection. In this paper, we aim to design a timetable from the perspective of the system cost of urban rail transit, which pays more attention to energy consumption. Since the costs of staff and rolling stock are considered in the model, the conducting of an energy-efficient timetable will be effective in the operation.
The timetable problem consists of optimizing the arrival time, dwell time, and departure time of trains at each station on the line under the time relationship and passenger demand constraints [4]. In particular, an energy-efficient timetable pays more attention to the running process in the section, which greatly affects the energy consumption and running time. There is a long history of research (Howlett et al. [5], Yang et al. [6], Canca et al. [4], Huang et al. [7]) focusing on energy-efficient train operation. It is always solving the optimal train control problem based on the Pontryagin Maximum Principle (PMP). Therefore, it is the primary step to determine the optimal driving strategy under the Pontryagin maximum principle. However, because of the non-linearity of the train control principle, it is hard to get a precise and optimal solution when modeling the timetable problem into a non-linear formulation. Although many studies employed the non-linear model with solutions based on the heuristic methods such as the intelligent Genetic Algorithm (GA) [8], it might fall into local optima or find an approximate optimum by heuristic methods. To fill this gap, we manage to make the following contributions: (1) We develop an approach with a mixed-integer linear programming formulation to optimize the supplement of running time in each track, arrival time, dwell time, and departure time at each platform, which aims at reducing the total energy consumption in the planning horizon.
A mixed-integer linear programming model is proposed by linearizing the objective and constraints, and it can be solved by the modern solver (Gurobi) easily.
(3) Then, from the perspective of the system cost, we develop an alternative objective based on the energy-efficient model, which gives an elaborate timetable for the conservation of system costs during the operation. (4) Moreover, real data from Beijing Changping Line is employed to validate the integrated approach, and the results are logically explained.
The remainder is organized as follows. The literature is reviewed in Section 2. Then, a mixed-integer non-linear programming model is developed in Section 3. By linearizing the objective and constraints, a mixed-integer-linear programming model is proposed in Section 4. From a system perspective, we give an alternative objective in Section 5. In Section 6, real data from Beijing urban rail transit are taken to validate the model. Finally, Section 7 gives the summaries of the results and suggestions for the future.

Literature Review
Energy-efficient train operations make an essential contribution to cutting down the energy consumption of the urban rail transit system [6]. There are three main types, including optimal train control, timetable optimization, and integrated optimization.
A collection of papers published after the 1990s was focused on optimal train control [9][10][11]. During the preliminary stage, Howlett (1994) [5] proposed a discrete model that considered the speed limits with fixed running time and used Kuhn-Tucker conditions to prove the optimal control strategy. Furthermore, Howlett and Cheng (1997) [12] solved the matter by applying continuously-varying gradients. Then, Khmelnitsky et al. (2000) [13] and Albrecht et al. (2013) [14] utilized independent arguments with different conditions and ensured that the optimal strategy was unique. Rodrigo et al. (2013) [15] discretized the speed profile and applied the Lagrange multipliers method to find semi-analytical solutions. Recently, aiming to find the optimal speed profiles and control regimes to cover the safety, punctuality, and energy efficiency of train operations, Wang et al. (2016) [16] studied the train trajectory optimization problem, which considered both the operational and signaling constraints. Paying attention to delay cases, with the two objectives of reducing delays and reducing energy consumption,   [17] proposed a multi-train trajectory optimization for single-track lines. However, the traditional optimal control method is complex and needs many physical parameters of the train that were tested by the train factory. Therefore, a data-driven optimization model emerges. Huang et al. (2019) [7] proposed a Data-Driven Optimization Model (DDOM) to describe the relationship between energy consumption and the discrete speed profile, then integrated two typical machine learning algorithms, Random Forest Regression (RFR) and Support Vector machine Regression (SVR), into a heuristic algorithm to solve the model.
Similarly, much research has been associated with the train timetable optimization problem in rail transportation systems (Shakibayifar et al., 2017 [18], Wang et al., 2017 [19], Gainanov et al., 2017 [20], Hassannayebi et al., 2018 [21]). These studies have enhanced the optimization theory and practice in the timetable problem. However, the timetable problem involves many sub-problem and still needs more research. Under the constraint of providing enough frequencies with acceptable services, Hu and Liu (2014) [22] proposed a headway-oriented timetable model to optimize the cost of a mass rapid transit system. Paying attention to the coordination between the last trains at the transfer stations of the urban rail transit networks, Kang et al. (2015) [23] proposed a method to reveal the relationships between the passenger transfer connection time and passenger transfer waiting time. To prevent passenger congestion and long waiting through train scheduling, Hassannayebi et al. (2016) [24] aimed to minimize the average waiting time of the passengers by using a path-indexed non-linear model to optimize the train timetable. Under the condition of emergencies, with the purpose of achieving the minimum delay time in the urban rail network, a binary mixed-integer programming model was presented by Huo et al. (2016) [25] based on the priority and train order entropy. To stress the level of transfer synchronization in the urban rail network, Guo et al. (2017) [26] presented a mixed-integer non-linear programming model with the objective of maximizing synchronization of the transfer process, which was solved through the combined heuristic of Particle Swarm Optimization and Simulated Annealing (PSO-SA). The energy consumption was described as a diminishing function of the section's running time. Canca and Zarzo (2017) [27] presented a mixed-integer non-linear model to optimize the timetable, which was solved using a sequential mixed-integer linear procedure. To handle the train departure orders, running time, dwell time, and overtaking locations flexibly, Sparing and Goverde (2017) [28] described an optimization method for designing a periodic timetable with the aim of finding the maximum stability within a heterogeneous railway network.
Integrated optimization in energy-efficient operations has recently become a hot issue, especially combining timetables with speed profile optimization (Zhou et [33] proposed a multi-objective evolutionary algorithm with the objectives of finding the optimal running time and reducing energy consumption. In addition, focusing on improving the utilization of train regenerative energy through the timetable optimization has also attracted more attention (Yang et al., 2016) [6]. Taking into account the dynamic train scheduling, Li and Lo (2014) [34] developed a control framework that simultaneously forecasted passenger demand, as well as optimized the timetable and speed profile. When the transportation service operator and passenger costs were simultaneously considered in a dense passenger Railway Rapid Transit (RRT) network, a Mixed-Integer Non-Linear Programming (MINLP) model was employed by Canca et al. (2016) [4] for optimizing the frequencies and determining the capacities. Time-dependent passenger demand in oversaturated conditions was also considered with the traffic capacity constraints. Then, in the operation scenario of two trains in multiple inter-stations, Huang et al. (2018) [2] integrated the operation strategies with optimizing the train timetable in pursuit of minimizing the energy consumption. Taking multi-phase speed limits into account, Yang et al. (2018) [35] employed an energy-efficient timetable model, which determined both the timetable and speed profile in each section by Particle Swarm Optimization (PSO) methods. Then, through the Taylor approximation, they transformed the train scheduling problem using a non-convex formation into a quadratic formation. Finally, an active set method was employed to search for the approximate optimal solution. Taking into account both the passenger waiting time and energy consumption under oversaturated conditions, Sun et al. (2019) [36] presented a bi-objective model and conducted a case study using a GA search in the Beijing subway system. By introducing a stop-skipping pattern before timetable optimization, Yang et al. (2019) [37] proposed a two-step optimization-based approach, which determined the specific skipped station based on the passenger smart-card data in the first-step, then modeled the energy-efficient timetable problem in convex quadratic programming, which was solved by the CPLEX. When the abnormal condition occurred (like delay perturbations), Yang et al. (2019) [38] employed an energy-efficient rescheduling approach to minimize the net energy consumption under the premise of reducing or eliminating the delay altogether.
Based on the review above, as we know, modeling the energy-efficient timetable problem by a special linearization technique to determine the timetable and train profile simultaneously has not been considered well. Furthermore, the complex model or using the heuristic algorithm to get an approximate solution should be improved. Therefore, this research is a supplement for relevant studies.

The Timetable Optimization Model: Mixed-Integer Non-Linear Programming Formulation
In this section, a timetable optimization model with a mixed-integer non-linear formulation is proposed to minimize the energy consumption in the planning horizon. The urban rail line is defined as a line network with a platform and track, which represent the node and link, respectively. The notations of the model are defined in the following.

Symbols
Parameters:  Decision variables: a p arrival time of the first train at platform p ∈ P, d p departure time of the first train at platform p ∈ P, dw p dwell time of the first train at platform p ∈ P, TB p boarding time at platform p ∈ P, TA p alighting time at platform p ∈ P, L t average train passenger load in track t, M t total train load in track t, E t energy consumption in track t, R t running time in track t, v t average velocity in track t, f frequency in the planning horizon, C circle time, h headway in the planning horizon.

Energy Consumption Function
A bi-directional line network was taken into account. There are N stations and 2N platforms, as shown in Figure  We can see that if station s ∈ S, then s ∈ P and 2N − s + 1 ∈ P. For each track t ∈ T, there are some attributes: length of track l t and maximum and minimum average velocity V max t and V min t (for the technology requirement and safety considerations of the rail line).
Up direction Down direction When the passenger demands and passenger-oriented objectives are considered to design a more energy-efficient timetable, both the departure and arrival time at each platform need to be determined. As the time is obtained, the trains will repeat the same running cycle with an interval of headway time (h), which is called the periodic timetable. The period usually lasts for 1-2 h. For instance, there is a cycle of 1 h in the evening peak from 17:30-18:30 in many metropolises. In the planning horizon, a p denotes the first arrival time at platform p ∈ P, while d p denotes the first departure time at platform p ∈ P, and v t represents the running speed at an average level on tracks t ∈ T. Therefore, the average running time of track t can be calculated using the formula R t = l t /v t . To pursue the optimal energy consumption, once the average speed is given in each track, the train will be operated under the optimal profile, as mentioned in Howlett et al. (2009) [39], with the Pontryagin maximum principle.
There are several parts of energy consumption during the running process: (1) Energy consumption that comes from the running resistance E rt . (2) Energy consumption that is from the air volume, which is absorbed by the train E it .
Energy consumption that comes from the aerodynamic resistance E at . (4) Energy consumption that offsets the kinetic energy, which is from the variation of the train on the tracks E kt .
Then, we can calculate the total energy consumption in the track t as follows: The total energy is a quadratic polynomial in v t , which is an increasing function of the passenger load and speed. In this study, the total energy with various passenger loads is assumed to be a direct ratio to its growth. Furthermore, we can obtain the basic energy consumption of an empty train from the real recorded data. That is to say, if we get the passenger load and train mass, we can calculate the total energy. As shown in Figure 3, there are real-world statistics for different sections (different colors and trend lines represent the energy consumption in different sections). It is easy to observe that the most energy-efficient timetable will be obtained under the minimum average speed (maximum allowed running time). However, this is a single train condition with no constraints. If the headway and passenger accommodation constraints are taken into account, then it would be more complicated. The minimum average speed can always be invalidated because of the pressure from the operational efficiency.
The energy consumed is a function of the average velocity (also the running time R t ) and train load M t (train mass M T and the total weight of passengers L t on the track t) in each track.
Furthermore, the energy consumption in a cycle is the sum of each track t ∈ T with the train running a whole circle on the bi-directional urban rail line. Therefore, it can be formulated as follows: Therefore, the total energy consumed during the planning horizon, E, including every train's energy in the cycle, should be computed as follows: We know that frequency f is associated with the train velocity since the faster a train runs, the higher the frequency is during the planning horizon, which also means higher energy consumption. However, if the frequency is too low, then the operational trains might fail to accommodate all passengers.

The Timetable Specification
OD: passenger demands in the form of the origin-destination (o-d) matrix: There are two conditions. In the upward direction, i < j, od i,j corresponds to the passenger volume from platforms i to j. Otherwise, for i > j, od i,j is the passenger volume from platform Thus, in the planning horizon (1 h in this paper), the boarding and alighting passengers can be calculated as follows: Firstly, the passengers boarding B p at platforms p can be computed as follows: For platforms N and 2N, there are no passengers boarding the train. Therefore: Similarly, the passengers alighting A p at platform p can be computed as follows: For Platforms 1 and N + 1, there are no passengers alighting from the train. Therefore: As for passenger volume in each track t, it can be computed through the accumulative process in linking platforms t to t + 1. For instance, in the upward direction (t < N), as shown in Figure 4, Therefore, the passenger volumes in both directions are: The largest passenger track can be found, and the maximum number of passengers in the tracks is: Q max is an essential parameter for accommodating the passenger demand. Furthermore, it is a constraint on transport capacity.
Assuming that the passengers arrive uniformly, the average passengers boarding and alighting the trains are applied in the model. A planning horizon of 1 h (3600 s) was taken into account in this study. Average boarding and alighting rates are (measured as passengers per second): The running time on track t begins at the departure time at platform t and ends at the arrival time at platform t + 1. Therefore: Furthermore, there are always velocity limitations in the running sections for safety and geographical reasons. Therefore, there are bounded constraints on the running time of sections: Similarly, the dwell time, dw p , begins at the arrival time at platform p and ends at the departure time of platform p. Therefore: In addition, there are constraints on the boarding and alighting of passengers during the dwelling process. Furthermore, a (r) denotes the empirical alighting, while b (r) denotes the boarding rates. Both are measured as seconds per passenger, which are determined by the width and number of train doors. Furthermore, the alighting time TA and boarding time TB at platform p can be computed as follows: The dwell time must be long enough to complete the alighting and boarding processes (as shown in Figure 5). Therefore, the dwell time should cover the alighting time and boarding time: dw p ≥ (TA p + TB P ), which are arranged as: There are always bounded values for headway. On the one hand, the minimum headway is determined by the signaling system and safety considerations. On the other hand, the maximum value is usually restricted to the service level or government tasks. Therefore: The departure frequency f , which represents the train quantity in the line during the planning horizon, is associated with the headway by: To ensure the frequency is an integer, h must be a divisor of 3600. In this case, as in Canca et al. (2016) [4], there are alternative headway values, for instance 120, 180, 240, 300, 360, 600, 720, 900, 1200, and 1800 s.
To avoid train pileups, the dwell time at each platform cannot be greater than the headway: Furthermore, to maintain the order of the passengers on the platforms for safety purposes, there is bounded dwell time: The cycle time C, which is the time consumed by a train in running a loop in the bi-directional line, can be computed as: where TZ denotes the turning back time. There is always a technical standard for the turn back time under the turn back facilities and signal system. Furthermore, the headway is associated with the fleet size F, which is a relationship in: Then, by replacing the headway h with f /3600, the fleet size can be calculated as: We can see that if cycle time C ≥ 3600 s, then F ≥ f . If cycle time C < 3600 s (there may be short lines for a short cycle), then F < f .
The train quantity is also limited by consideration of the fixed investment. Therefore, the fleet size cannot be greater than the maximum train quantity: For a fixed train unit, there is a corresponding capacity Train cap due to the consideration of the limited room and passenger comfort. To guarantee that all passenger demands (all OD demands) are satisfied with the maximum capacity, frequency f should satisfy the constraint as follows: The constraints (23) can be another form by replacing the f by relation f × h = 3600, then: Therefore, if the train capacity is fixed after determining the train units, then the carrying quantity in one hour is Train cap /h min . Additionally, with the increase in Q max , the headway should be shorter.
Finally, if the average passenger weight is m 0 , the average train passenger load L t on track t for a train in the planning horizon can be calculated as:

First MINLP Formulation
Therefore, the model is a Mixed Integer Non-Linear (MINLP) problem. Moreover, non-linearities include the objective function E and the constraints ( †) and ( ‡). Since we are in pursuit of the exact solution and a global optimum, we will linearize them in the following section.

The Timetable Optimization Model: Mixed Integer Linear Programming Formulation
To find the exact solutions more easily, we linearize the objective function and non-linear constraints in this section and develop a timetable optimization model with a mixed-integer linear programming formulation. The linearizing process is concluded in the following.

Linearizing the Objective Function
In real operations, the profiles in a section are not arbitrary, but a set of alternative train speed profiles. There are different operation levels g, which correspond to different running times R g . When the profile level is changed, the running time and energy consumption will be changed. The basic energy consumption of an empty train on track t at level g is E 0 t,g . For the sake of clarity, as shown in Figure 6, there are three levels to the operational process in the section, which is from Ming Tomb to Changpingxishankou of the Changping Line in Beijing urban rail transit. Thus, there are three profiles with three different running times and corresponding energy consumption. When the passengers are loaded on the train, then the total weight is M T + L t . Thus, the energy can be approximately equal to (1 + L t M T )E 0 t,g . Consider that there is a set of operation levels on track t, G t . Introduce a set of binary variables θ g,t . If the level g is selected, then the θ g,t = 1; else θ g,t = 0. For a running section, there is only one level g that can be selected; therefore: Therefore, the energy consumption E t (R t , L t ) can be replaced by: Define a new set of variables ρ g,t so that ρ g,t = (1 + L t M T ) if θ g,t = 1. Then, a set of constraints should be introduced: ρ g,t ≤ θ g,t Λ, ∀t ∈ T, g ∈ G t (31) ρ g,t ≥ 0, ∀t ∈ T, g ∈ G t (32) where Λ is the upper bound; for instance, (1 + L max M T )E 0 max . Therefore, the objective function E t (R t , L t ) can be replaced by: We introduced a series of binary variables δ m (m = 1, 2..., |Ω|). Since the passenger demand is strong in the peak hour, Ω = {30, 20, 15, 12, 10, 6} and H = {120, 180, 240, 300, 360, 600} are in our consideration. The binary variables δ m and the elements h m can be utilized to represent the headway as the sum of their product as follows: Similarly, the frequency can be expressed as: Therefore, the sum of the binary variables δ m should be equal to one: By replacing with the new variables, the objective can be expressed in the following formation: In the above formulas, replace (ρ g,t × δ m ) with the new variable λ g,t,m , which must take the value ρ g,t while δ m = 1. This can be achieved through sets of constraints: where ∆ is the upper bound; for instance, (1 + L max M T )E 0 max f max . Finally, the linear formulation of the objective is shown as follows:

Linearizing the Non-Linear Constraints
The constraints † and ‡) are non-linear. After introducing the binary variables δ m and sets of H and Ω above, the constraint † is linearized. Then, the constraint ‡ can be rewritten in: Then, a set of variables α m is introduced. If δ m = 1, then α m = F. Therefore, it can be achieved by sets of constraints: where Ψ is the upper bound; for example, Ψ = F max × h max . Then, the constraints ‡ can be rewritten as follows:

An Alternative Objective Function from the Perspective of System Cost
The aim of designing a timetable is optimizing energy consumption, which is pretty significant for environmental protection. However, it seems not profitable for urban rail system operators. The reason is that there are other costs like driver salaries, train unit acquirement, maintenance cost, etc., during the operations. Among these costs, two main ones come from the personnel expenses and those generated from train unit acquisition. Therefore, to minimize the total cost from an operational system, we employed an alternative objective function, which takes into consideration the train unit and operation cost. There are three types of value coefficients as follows: C e : The average electricity price of 1 kWh. C a : The hourly cost, which comes from train unit acquisition and maintenance during the operations. C s : The average hourly salary corresponding to the driver. Therefore, the alternative objective function can be: We can see that the first term of the sum is the energy cost, and the latter corresponds to the cost of investments in the train unit and personnel.

A Case Study: Based on the Beijing Changping Line
To give a detailed illustration of the proposed model, we considered the bi-directional line (Changping Line) of Beijing urban rail transit. The operational section of Changping Line is a two-way line running from the Xi'erqi station to the Changpingxishankou station, with operational mileage of 31.9 km and a total of 12 stations opened (as illustrated in Figure 7). Since the Changping Line is located in the suburbs of Beijing, it is not overcrowded in the peak hours. Information about the stations, the length of open tracks, ranges of the dwell time, and profile attributes are given in Table 1. The profile energy consumption was derived from the real experiment data of the empty train. The dwell time limits were in consideration of the passenger waiting time and maintaining platform order.
The operational parameters are listed in Table 2. There are 22 trains belonging to the operation branch of Changping line. The net weight of each train is about 205 ton, while the train seating capacity is1460 people. The peak hour coefficient is about 1.2 based on a survey of the actual ride. Therefore, the train capacity is 1.2 * 1460 = 1760 people. The average passenger weight is given as 65 kg. The alighting and boarding rates are empirical, which refer to Canca et al. (2017). The turning back time for the Beijing urban rail transit networks is a technical standard that should be no more than 5min. The cost coefficients are obtained from the operations department. Table 3 shows the passenger matrix for the peak hour of the planning horizon 7:00-8:00, which was extracted from the AFC (Automatic Fare Collection) data of Changping Line. The passenger demand (boarding and alighting numbers) and load profile at each station are shown in Figure 8. We can see that the maximum section volume (Q max ) was from Beishaowa-Changpingdongguan, with about 22,156 passengers. We can see that there was a slight imbalance in passenger demands between the two directions, which was a result of more passengers commuting from a residential zone to work places or schools.
Case 1: The objective function includes minimizing the energy consumption only. We obtained the minimum energy consumption (9413.3 kWh) and timetable under all constraints. To show the energy savings, we calculated the maximum energy consumption (14,458.5 kWh, which was under the fastest profile. Therefore, 14,458.5 − 9413.3 = 5045.2 kWh were saved, which accounts for about a 35% reduction. The results of arrival time, departure time, and dwell time of the first train in the planning horizon [0, 3600] s are listed in Table 4. We can see that the dwell time was relatively even, except for three platforms, which is in accordance with the passenger flow. Therefore, in most of the stations, the passenger flow was not overly crowded. For clarity, all departures in both directions are shown in Figure 9. It illustrates that the distributions of the departure time in the up and down direction was symmetrical. Furthermore, the arrival time in the up direction was earlier than that in the down direction, which may be caused by the imbalance in passenger demands between the two directions; because the passenger flow in the up direction was higher than that in the down direction.   The running time and energy consumption of each track are shown in Figures 10 and 11. Compared with the max running time of each track, the actual running time was equal to those of most tracks. Correspondingly, the actual energy consumption levels were equal to the minimum of most tracks, which conforms to the pursuit of minimum energy consumption. Furthermore, we calculated the travel time of the whole odpairs as shown in Figure 12. There were two surfaces of the up and down direction. The surfaces were relatively smooth, which indicates that the running time and departure time were well distributed. From the intersection angles with their horizontal plane, we can see that the travel time in up direction (from Station 1-12) increased slightly more quickly than in the down direction (from Station 12-1). Overall, the distributions of travel time were balanced in the whole line. was 6.6% compared with the maximum system cost 55,880.9 (RMB). The frequency and headway were the same. Therefore, there was little difference in the waiting time of the passenger.  Running time and energy consumption with two different objectives are given in Table 6. Thus, we see that the running time under the alternative objective in each track tended to be shorter than that under the energy objective. Therefore, from the viewpoint of the system cost, the train tended to run faster and made full use of the train fleet size. For clarity, as shown in Figures 13 and 14, compared with the maximum time and minimum energy consumption level, the actual running time under the alternative objective deviated from the maximum time, which shows high efficiency in the operation. In summary, the results of the alternative objective were a compromise between the maximum energy consumption and energy objective.    Table 6. Running time and energy consumption with two different objectives.

Items Running Time (s) Energy Consumption (kWh)
As shown in Figure 15, the departure time of the train under the alternative objective was earlier. In addition, the high operational efficiency may cause extra energy consumption, but a reduction in total cost. Therefore, there is a trade-off between energy conservation and operational efficiency. The solutions under the balanced condition between energy conservation and operational efficiency were found by the model. Furthermore, whether the view was oriented toward energy conservation or system cost, the departure times in the up and down directions were in equilibrium.

Conclusions
Urban inhabitants relying too much on urban rail transit due to its advantages of safety, punctuality, high energy, and large capacity, there has been a quick rise in energy consumption. For the sake of pursuing cost savings and environmental protection, there is an urgent demand for energy conservation.
To address the problem of increasing energy consumption in the urban rail transit system, we aimed to design an energy-efficient timetable in this paper. First, we proposed a timetable optimization model in the mixed-integer non-linear formulation. Then, the objective and constraints were linearized. The mixed-integer linear programming model was employed, which could be solved by the commercial solver Gurobi. Moreover, from the perspective of system cost, we presented an alternative objective. To validate the model, real data from the Changping Line of Beijing urban rail transit were applied as a case study. Subsequently, the results showed that the designed timetable could achieve about a 35% energy reduction compared with the maximum energy consumption and 6.6% cost savings compared with the maximum system cost. We find that there was a trade-off between operational efficiency and energy consumption. The energy objective was prone to lower energy consumption, while the alternative objective tended toward high operational efficiency with fast running time. Furthermore, the direction with higher passenger flow also tended to run faster.
The modeling framework still had some limitations. For instance, the passenger flow was roughly considered. A dynamic passenger flow will be more practical. In summary, there are several future research directions. We could consider passenger behavior and passenger costs, such as the travel times of passengers. There is feedback about passenger flow after a timetable change. Therefore, a dynamic equilibrium between the timetable and passengers should be found. Moreover, the regenerative energy produced during the braking process could be utilized in the same substation, which could reduce both the energy and total costs, so it should be modeled. Finally, there is not always a single line, but an urban rail transit network in a city. The energy conservation should be extended to the urban rail transit networks for the entire operational system. There are several future research directions. Finally, a scientific algorithm should be proposed when this problem is for a large-scale condition.
Author Contributions: All the authors contributed to the work of this paper. H.L. and K.H. designed and developed the main parts of the research work, including the design of the model, writing the original draft, and programming. Y.Z. and J.W. contributed to the writing, review and editing. X.Y. was actively involved in verifying the work and finalizing the manuscript. All authors worked together on this paper, and all authors read and approved the final manuscript.