Preprint
Article

Energy-efficient train driving based on optimal control theory

Altmetrics

Downloads

155

Views

72

Comments

0

A peer-reviewed article of this preprint also exists.

This version is not peer-reviewed

Submitted:

10 August 2023

Posted:

10 August 2023

You are already at the latest version

Alerts
Abstract
Efficient train driving plays a vital role in reducing the overall energy consumption in the railway sector. An energy minimising control strategy can be computed using the framework given by optimal control theory, in particular the Pontryagin maximum principle. Our optimisation approach is based on an algorithm presented by Khmelnitsky that considers electric trains equipped with regenerative braking. A derivation of Khmelnitsky’s theory from a more general formulation of the maximum principle is given in this article, and a complete list of switching cases between different driving regimes has been included that is essential for practical application. A number of numerical examples are added to visualise the various switching cases. Energy consumption data from real-life operation of passenger trains are compared to the calculated energy minimum. In the presented study, the optimised strategy was able to save 37 percent of the average energy demand of the train in operation.
Keywords: 
Subject: Engineering  -   Transportation Science and Technology

1. Introduction

Energy-efficient driving strategies have shown to substantially lower the total energy consumption of railway trains. With an increased demand for energy saving technologies, numerous research studies have been devoted to this subject in recent years, showing that in many cases an average energy saving of 15 to 30 % could be achieved by assisted driving [2,3,4,5,6,7,8,9]. A comprehensive review on energy-efficient train control, commonly abbreviated as EETC, is given by Scheepmaker, Goverde, and Kroon [10]. The majority of the mathematical models that have been developed to calculate energy minimised driving are based on Pontryagin’s maximum principle (PMP), a fundamental theorem of optimal control theory. The maximum principle was formulated in 1961 by Pontryagin, Boltyanskii, Gamkrelidze, and Mishchenko [11,12]. The first application of Pontryagin’s maximum principle to the operation of trains dates back to Japanese research in 1968 (Ichikawa [13]). The model of Ichikawa was restricted to a flat track and linear train resistance depending on speed. Another early contribution is due to Strobel, Horn, and Kosemund [2]. They considered non-constant altitude and found that the optimal solution is always in one of the drive regimes full traction, partial traction with constant speed, coasting, partial braking with constant speed, and full braking. Based on their model, they set up a first driver advisory system for the Berlin S-Bahn (suburban trains).
Since those pioneering research on this topic, models have gained complexity. Important developments have been: (1) the application to variable altitude and variable speed limit, (2) the inclusion of regenerative braking in electric locomotives or motor coaches, (3) improved modelling for train resistance, tractive, and brake force depending on speed, (4) the energy minimisation with non-zero speed boundary conditions, (5) advanced modelling of efficiency as dependent on speed, tractive, and brake force, (6) multiple-train problems.
Regarding the second aspect, some (especially the earlier) authors are dealing only with mechanical braking of trains. When regenerative braking is considered, most authors base their energy-efficient driving strategy on the exclusive use of the regenerative brake. A few recent contributions cover the combined application of mechanical and regenerative brake in their optimisation strategy.
We proceed with some important developments for purely mechanical braking. Starting from the PhD thesis of Milroy [14] in 1980, a reaserch team at the University of South Australia introduced the system Metromiser (Howlett and Pudney [3]). The system covered both timetable planning and driver advisory. Metromiser was restricted to constant track gradients during both coasting and braking. The system was first applied on trains in Adelaide (Australia) in 1984, followed by Toronto (Canada), Melbourne (Australia), and Brisbane (Australia). Howlett, Milroy, and Pudney [15] were the first to assume variable speed limit, but only for the special case of discrete traction control. In 2003, Liu and Golovitcher [16] presented a method covering variable altitude, variable speed limit, continuous traction control, and quadratic train resistance dependence on speed. Except being bound to purely mechanical braking, this algorithm is therefore quite general with respect to the other assumptions. The article contains a complete description of switching between the possible drive regimes. In 2013, Albrecht, Howlett, Pudney, and Vu [17] introduced the driver advisory system Energymiser that has been applied by French railway SNCF on TGV high speed trains. They proved the uniqueness of the energy optimum being derived. A multi-objective algorithm combining energy minimisation and punctuality of trains was developed by Aradi, Becsi, and Gaspar [7].
Regenerative braking was first considered in 1985 by Asnis, Dmitruk, and Osmolovskii [18] in the context of train energy minimisation. They assumed constant altitude and speed limit in their model. Franke, Terwiesch, and Meyer [4] introduced an energy optimisation algorithm for variable altitude, variable speed limit, and purely regenerative braking that was restricted to piecewise constant traction and brake force. They suggested a discrete dynamic programming (DDP) method for solution. In 2000, Khmelnitsky [1] developed a model that covered variable altitude, variable speed limit, continuous traction control, arbitrary maximum tractive and brake force depending on speed, quadratic train resistance with respect to speed, and exclusively regenerative braking. Regarding this features, it is a broadly applicable approach with little restrictions. Due to this advantages, the algorithm of Khmelnitsky has been chosen as the basis for our research. In a recent article by Ying et al. [19], the EETC problem solved by Khmelnitsky has been revisited, with a special focus on solutions touching the speed limit. The paper of Ying et al. [19] can be especially recommended for its comprehensive illustration of the large variety of cases that can arise during construction of the optimum solution. The EETC problem with prescribed non-zero speed boundary conditions was solved by Ying et al. [20].
There are relatively few research contributions dealing with the EETC problem using combined mechanical and regenerative braking. Baranov, Meleshin, and Chin’ [21] were the first considering this topic, but an algorithm to solve the problem was left as an open question in their article. Lu et al. [22] studied combined mechanical and regenerative braking, but excluded both variable speed limit and variable track gradients. Zhou et al. [23] studied synchronisation of accelerating and braking trains, considering both kinds of braking, but again did not include variable speed limit. Fernández-Rodríguez et al. [24] combined both brake systems in a multi-objective optimisation method, but did not derive a rigorous energy minimum. In the paper of Scheepmaker and Goverde [25], the optimisation problem with combined mechanical and regenerative braking was solved for the first time under general assumptions, i.e., with speed dependent tractive and brake effort, variable speed limit, and variable altitude. Scheepmaker and Goverde have used the Gauss-Radau pseudospectral method implemented in GPOPS [26] to solve the optimal control problem.
Several authors have considered advanced efficiency modelling, in particular efficiency depending on speed, tractive, and brake force. This assumption is more realistic but substantially complicates the approach. Most notably, the optimal trajectory for this problem is no more restricted to the driving regimes of full traction, constant speed, coasting, and full brake. Therefore, other methods than PMP techniques are used in this case, in particular track length discretisation combined with nonlinear programming. Research on advanced efficiency modelling was carried out, e.g., by Franke, Terwiesch, and Meyer [4], Ghaviha et al. [27], Kouzoupis et al. [28], and Feng, Huang, and Lu [29].
Multiple train problems have often been studied with a focus on timetable design in order to synchronise accelerating and braking trains for best distribution of regenerative energy. See, for example, Zhou et al. [23] for research on this topic. More references can be found in Scheepmaker and Goverde [25]. Szkopiński and Kochan [30] have studied energy-efficient train driving when approaching a train on the track ahead.
Apart from the many research contributions based on an application of Pontryagin’s maximum principle, a variety of other methods exist to solve the EETC problem. In particular, so-called direct methods have recently gained attention. They are characterised by the fact that the problem is first discretised, usually by dividing the track length into smaller intervals. The resulting nonlinear programming (NLP) problem is then solved by nonlinear programming methods as, for example, the pseudospectral one. Wang et al. [31] were the first to apply the pseudospectral method to an EETC problem. Scheepmaker and Goverde [25] also utilised a pseudospectral approach. Direct methods have the advantage of being very flexible with respect to the problem under consideration. In addition, optimal control solvers like DIDO (Direct and Indirect Dynamic Optimization) [32,33], PSOPT [34], or GPOPS (General Pseudospectral Optimal Control Software) [26] are readily available. The implementation is much easier than a solution of the problem along the lines of Pontryagin’s maximum principle. On the other hand, the pseudospectral method often shows an inaccurate oscillatory behaviour of solutions, and it tends to be time consuming. [28] Recently, Kouzoupis et al. [28] were able to reduce the computing time of a direct method using multiple shooting.
Based on the above-mentioned approach of Khmelnitsky [1], a programme called opTop (optimum train operation) has been developed at Fraunhofer Institute for Factory Operation and Automation (IFF), Magdeburg, Germany. The code is written in MATLAB and currently features energy minimisation with exclusively regenerative braking, except for fastest train motion in the case of tight timetables where additional mechanical braking is considered. The choice of Khmelnitsky’s method was mainly motivated by superior accuracy at an acceptable computing time. An extension of the code to fully include mechanical braking into the optimisation is subject to further development, as well as non-zero speed boundary conditions and smooth switching from fastest to energy-efficient driving in case of a train delay. The code has been tested offline in numerous cases based on real railway tracks and timetables, and computing time has been substantially reduced by code optimisation. A couple of tests with driver advisory in real-life operation have shown energy savings of about 20 % compared to an average of unassisted runs.
The present article was driven by two main motivations. The first one was to show a direct derivation of Khmelnitsky’s theory from a more general formulation of Pontryagin’s maximum principle given by Hartl, Sethi, and Vickson [35]. By seeing Khmelnitsky’s theory in this general framework, extensions to other conditions could be elaborated. A second motivation was to provide a comprehensive illustration of the behaviour of the method of Khmelnitsky using a couple of numerical examples. We have felt that this very useful method would strongly benefit from some illustrative examples that show the behaviour of the trajectory field including kink points as well as the large variety of switching cases from one driving regime to another. We have, therefore, put a strong focus on the examples in Section 2.4 and Section 2.5, in which most of the possible switching cases can be studied in detail.
The article is structured in the following way: In Section 2.1, the model of train motion is introduced. Section 2.2 is devoted to the Pontryagin maximum principle that provides a couple of necessary conditions the energy minimising solution has to fulfil. In Section 2.3, the maximum principle is applied to the specific minimum energy problem for train motion. This leads to the observation that only four driving regimes full traction, constant speed, coasting, and full regenerative braking are feasible for an energy minimising motion. The regime constant speed can be driven only at two specific velocities or at speed limit. Section 2.4 and Section 2.5 describe the algorithm of Khmelnitsky [1] for the construction of the minimum energy solution. A complete list of switching cases is given and explained in Section 2.4, and four numerical examples are introduced to illustrate those cases. Some directions are given that lead to a substantial reduction of computing time. Finally, in Section 3, the energy optimum is compared to energy consumption of various train runs on a specific railway line in real-life operation. As a result, the optimal strategy consumed only 63 % of the average energy demand of the unassisted runs in real operation.

2. Materials and methods

2.1. Model of train motion

The equation of motion of a train with an electric engine and both regenerative and mechanical brake can be given by
F tr F br F mbr F air F roll F sl = c rot m a ,
where F tr is tractive force, F br is the force applied by the regenerative brake, F mbr is the force of the mechanical brake, F air is air drag, F roll is rolling friction, F sl is the downhill slope force, c rot is a coefficient representing rotating masses, m is the mass and a the acceleration of the train. This is Newton’s second law of motion with the extension that rotating masses as, for example, the wheels of the train, are accounted for by a factor c rot . A model of this type is commonly be used for train motion, see for example the review article of Scheepmaker, Goverde, and Kroon [10] or the textbook of Ihme [36]. According to Ihme [36], c rot = 1.06 1.11 for passenger trains, depending on their length, where longer trains will generally have smaller values of c rot .
The tractive force F tr is limited by both engine power and rail friction (adhesion). According to Fassbinder [37],
F tr F tr , max = min ( P tr , mech / v , μ ad g m loc ) ,
where P tr , mech is the mechanical engine power used for traction, v is the velocity of the train, μ ad is the adhesion coefficient, g is gravity, and m loc the mass of the locomotive. The adhesion coefficient μ ad has been obtained experimentally in 1943 by Curtius and Kniffler [38], see, e.g., Schlecht [39], leading to the empirical relation
μ ad = 7.5 v k m / h + 44 + 0.161 .
The adhesion coefficient μ ad attains a maximum value of 0.331 when v = 0 .
The regenerative brake force F br is, as the tractive force, restricted by engine power and rail adhesion. However, it must be further limited to avoid a derailing of coaches behind the braking locomotive. Thus, for the regenerative brake force there holds
F br F br , max = min ( P br , mech / v , μ ad g m loc , F br , lim ) ,
where P br , mech is the mechanical engine power for regenerative braking, and F br , lim is the additional limit to avoid derailing. In general, P tr , mech = P br , mech will hold. According to [40], the limit force F br , lim has been recently enlarged in Germany from 150 k N to 240 k N , a value that is also considered in Scheepmaker and Goverde [25]. The maximum forces F tr , max and F br , max are shown in Figure 1.
Remark. In practical operation, regenerative braking of trains is not possible when the speed is too small, meaning that for a final halt, the mechanical brake always has to be applied. This has been pointed out by Scheepmaker and Goverde [25], who have estimated a minimum speed of 8 k m / h for application of regenerative braking based on data from Netherlands Railways. However, since kinetic energy below 8 k m / h is relatively small, we have neglected this consideration in our model.
When a train brakes, using the regenerative brake is clearly advantegeous with respect to saving energy. However, the force of the regenerative brake is limited and, especially for long freight trains, much weaker than the force of the mechanical brake. This is due to the fact that regenerative braking applies only to the wheels of the locomotive, while the mechanical brake acts on every wheel of the train. If the time required by timetable is too short for a certain distance, purely regenerative braking might not be sufficient. But, as Scheepmaker and Goverde [25] have proven by means of optimal control theory, the mechanical brake is always ‘second choice’ in an energy-minimising solution. This means that mechanical braking is only applied when regenerative braking is operating at maximum force, i.e., when F br = F br , max holds. In the optimisation method presented here, the mechanical brake is not incorporated into the theory and will only be considered in the calculation of the fastest possible motion. This will be explained in more detail in the following sections.
The mechanical work applied for traction is given by
W tr = s 1 s 2 F tr d s = t 1 t 2 P tr , mech d t ,
where s is track length and t is time. The electric energy required for traction is
E tr = W tr / η tr .
The efficiency η tr of electric locomotives is usually given in a range of 83 to 87 percent [25,37,41]. Within this paper, we will use an intermediate value of η tr = 0.85 . Likewise, the mechanical work applied for regenerative braking is
W br = s 1 s 2 F br d s = t 1 t 2 P br , mech d t ,
and the electric energy returned will be
E br = η br W br .
We use a braking efficiency η br = η tr = 0.85 according to Fassbinder [37]. The difference E net = E tr E br is called the net energy, which should be minimised.
In real operation, additional energy E add will be required that is not directly related to traction or brake, including energy for air commpression, air conditioning, lighting and so on. However, since this additional energy is mainly a function of time, it can not be reduced by a driving strategy when the total time is constant as given by the train’s timetable. Therefore, the additional energy does not enter the net energy minimisation problem and can be excluded here.
For air drag, Ihme [36] recommends the so called Hanover Formula that goes back to Voß, Gackenholz and Wiebels [42]. It takes the form
F air = 1 2 ρ air A ref c w v 2 ,
with air density ρ air = 1.25 k g / m 3 . In this formula, A ref is not the cross-sectional area of the train, but a reference area of 10 m 2 . The air drag coefficient c w is calculated according to
c w = c w , loc + c w , first + ( n 2 ) c w , middle + c w , last ,
n being the number of coaches. Ihme [36] gives the values c w , loc = 0.26 , c w , first = 0.13 , c w , middle = 0.10 , c w , last = 0.23 as appropriate for Intercity coaches.
The rolling friction can be modelled by F roll = c roll m g with c roll = 0.0015 according to Ihme [36]. The downhill slope force is given by
F sl = m g d z d s ,
where z is the altitude of the track.

2.2. The maximum principle

The problem of minimising the net energy of a scheduled train can be formulated as an optimal control problem, as it has been proposed by Khmelnitsky [1] who also presented an algorithm to obtain the unique net energy minimum. The algorithm of Khmelnitsky is essentially based on the maximum principle of optimal control theory that has been developed by Pontryagin, Boltyanskii, Gamkrelidze, and Mishchenko [11,12]. A variety of formulations of the maximum principle can be found in Hartl, Sethi, and Vickson [35]. While some of the results of the maximum principle for the particular problem of a minimum energy train ride are already given in Khmelnitsky’s paper [1], a direct derivation from the more general formulation of the maximum principle as presented in Hartl, Sethi and Vickson [35] is not included in Khmelnitsky. We are, therefore, going to show this derivation here. In Section 2.2, the maximum principle will be presented, based on the ‘Informal Theorem 4.1’ in [35]. The subsequent Section 2.3 contains the application to the train problem. Equipped with this general framework, extensions of the theory (for example the inclusion of the mechanical brake) are possible.
We consider the following problem: Let s be the track length coordinate between two stations at s = 0 and s = S . The motion of the train is described by the time t ( s ) it takes for the train to move from the first station to position s. We assume t ( 0 ) = 0 and a fixed duration t ( S ) = T given by the timetable. The velocity v of the train is restricted by a piecewise constant function v max ( s ) that accounts for speed limits on the track. Tractive and regenerative brake force, F tr and F br , are limited according to (2) and (4), respectively. The optimisation problem is to find a motion, i.e., a function t ( s ) , that minimises the net energy E net .
We introduce some notation, closely following Khmelnitsky [1]. Let E kin = c rot m v 2 / 2 be the kinetic energy of the train. Dividing by c rot m leads to a specific kinetic energy
K = E kin c rot m = v 2 2
limited by speed restrictions: K ( s ) K max ( s ) = ( v max ( s ) ) 2 / 2 . Likewise, a specific potential energy is defined by
P = E pot c rot m = g z c rot .
A total specific mechanical energy is then given by E = K + P . In the same way, specific forces are defined as
u tr = F tr c rot m ,
u br = F br c rot m ,
g tr = F tr , max c rot m ,
g br = F br , max c rot m ,
w = F air + F roll c rot m .
A total recuperation efficiency is given by α = η tr η br . The aim of minimising the net energy E net is equal to maximising
J = η tr c rot m E net = 0 S u tr ( s ) α u br ( s ) d s .
With this notation, the optimisation problem can be written in the following canonical form:
  • state vector x = ( E , t )
  • control vector u = ( u tr , u br )
  • state differential equations
    d E d s = u tr u br w
    d t d s = 1 2 K
  • state boundary conditions
    E ( 0 ) = P ( 0 )
    t ( 0 ) = 0
    E ( S ) = P ( S )
    t ( S ) = T
  • maximum
    J = 0 S u tr ( s ) α u br ( s ) d s max
  • state constraint E K max + P
  • control constraints 0 u tr g tr and 0 u br g br
Remark. Further constraints could be added here, especially 0 t T and K 0 . However, K > 0 will always hold for the energetic optimum solution inside the intervall ( 0 , S ) . The condition 0 t T follows from (21), (23), and ( 23 ) . Therefore, conditions 0 t T and K 0 need not be set explicitly as a constraint. Simultaneous traction and brake, i.e. u tr ( s ) > 0 and u br ( s ) > 0 for the same s, is also not explicitly prevented by a constraint since such a solution will surely not be energetically optimal.
Remark. In equation (20), we have excluded mechanical braking. This means that, in accordance with Khmelnitsky [1], we formulate the optimisation problem only for purely regenerative braking. An extension of the presented theory to the case with simultaneous application of both regenerative and mechanical brake is possible in a straightforward manner, and is, in a slightly different formulation, given by Scheepmaker and Goverde [25].
Introducing the notation
f ( x , u , s ) = u tr u br w 1 / 2 K ,
F ( u ) = u tr + α u br ,
g ( x , u , s ) = u tr u br g tr u tr g br u br ,
h ( x , s ) = K max K ,
we have an optimisation problem of the form given in Hartl, Sethi and Vickson [35]:
d x d s = f ( x , u , s ) ,
J = 0 S F ( u ( s ) ) d s max ,
g ( x , u , s ) 0 ,
h ( x , s ) 0 .
In Hartl, Sethi, and Vickson [35], the problem is formulated with a time t being the independent variable instead of s. This is due to the fact that many practical optimisation problems are formulated in a time-dependent way. In our case, however, the position-dependent formulation has some advantages, especially since also the speed restriction depends on position, not on time.
The function
H = F + λ T f
with Lagrange multipliers λ ( s ) is called the Hamiltonian of problem (31)-(34). (Here and in the following, ( · ) T means the transposed of a column vector, i.e. λ T f is the scalar product of vectors λ and f .) The Lagrange multipliers λ ( s ) are also called the costates of the problem. In order to agree with the notation of Khmelnitsky [1], we denote the costates according to
λ ( s ) = ψ ( s ) ψ T ( s ) .
Furthermore, the Lagrangian is introduced by
L = H + μ T g + ν h ,
where μ ( s ) and ν ( s ) are also called Lagrange multipliers. The Lagrange multipliers λ ( s ) , μ ( s ) and ν ( s ) are continuous functions in s, except for positions s i where h ( x , s i ) = 0 .
The ‘Informal Theorem 4.1’ of Hartl, Sethi, and Vickson [35] states the following necessary conditions for ( x , u ) to be a solution of the optimisation problem (31)-(34):
  • The control vector u ( s ) maximises the Hamiltonian H pointwise for any s [ 0 , S ] .
     
  • Their holds
    L u = H u + μ T g u = 0 .
  • The costate equation
    d λ d s = L x
    is satisfied.
  • Let μ i be the components of vector μ , and g i be the components of vector g . There holds
    μ i ( s ) 0 and μ i ( s ) g i ( s ) = 0
    for all i and all s [ 0 , S ] . This is called complementary slackness in Khmelnitsky [1].
  • Another complementary slackness condition,
    ν ( s ) 0 and ν ( s ) h ( s ) = 0 ,
    holds for all s [ 0 , S ] .
  • The following jump condition is satisfied: At every point s i where λ is discontinuous, there exists a number η i with
    λ ( s i + ) λ ( s i ) = η i h x ,
    H ( s i + ) H ( s i ) = η i h s ,
    η i 0 and η i h ( s i ) = 0 .
    Here, the argument s i corresponds to the left-hand limit, and s i + to the right-hand limit at s i . The vector equation (43) is meant component-wise.
The Hamiltonian maximisation condition (38) is called the Pontryagin Maximum Principle, and equations (39)–(42) are referred to as the Karush-Kuhn-Tucker (KKT) conditions [43,44]. In Hartl, Sethi, and Vickson [35], the maximum principle is presented also for the case of multiple state constraints, i.e., when h is extended to a vector. The maximum principle does not, as we shall see, tell the solution to the energy minimisation problem directly, but it provides essential information such that a construction of the solution becomes possible in an iterative trial-and-error process.

2.3. Application of the maximum principle to the energy minimisation problem

From the costate equation (40) it follows for the second component of λ
d ψ T d s = L t = 0 ,
since L does not depend explicitly on t. Equation (43) gives, again for the second vector component,
ψ T ( s i + ) ψ T ( s i ) = η i h t = 0 .
Thus, the Lagrange multiplier ψ T is a constant. Equation (39) results in
L u tr = 1 + ψ + μ 1 μ 3 = 0 ,
L u br = α ψ + μ 2 μ 4 = 0 .
From the costate equation (40), first component, we have
d ψ d s = L E = ψ d w d K + ψ T · ( 2 K ) 3 / 2 μ 3 d g tr d K μ 4 d g br d K + ν .
When travelling below the speed limit, i.e. h > 0 , then ν = 0 holds due to (42). (In Khmelnitsky’s notation, μ 3 and μ 4 are called a tr and a b , respectively.) We now apply a case distinction to the multiplier ψ , and evaluate u tr , u br , and the μ i using the conditions (38), (41), (48), (49), and the continuity of the multipliers. The results are given in Table 1.
From the values given in Table 1 it follows that, in the interior of the state domain, μ 3 and μ 4 can be directly expressed in terms of ψ :
μ 3 ( s ) = ψ ( s ) 1 i f ψ ( s ) 0 0 i f ψ ( s ) < 0
μ 4 ( s ) = α ψ ( s ) i f ψ ( s ) α 0 i f ψ ( s ) > α
The cases indicated in Table 1 are the only driving modi that are possible for the optimum solution of the problem. This means that the interval [ 0 , S ] can be completely segmented into subintervals S i = ( s i 1 , s i ) with 0 = s 0 < s 1 < < s n = S , where every subinterval S i corresponds to one of the drive modi full traction, partial traction, coasting, partial regenerative brake, and full regenerative brake. The cases full traction, coasting, and full regenerative brake are regular in the sense that the control variables u tr and u br are both defined, and therefore the equation of motion is completely given by the underlying physics. On the contrary, the modi partial traction and partial brake are singular, meaning that here only one control variable is defined. Therefore, they need some additional consideration. We will distinguish the following four cases:
Case 1: partial traction below speed limit. Let S i be an interval of partial traction in the interior of the state domain, i.e., an interval where ψ ( s ) 1 and v < v max hold. It follows μ 3 ( s ) μ 4 ( s ) 0 on this interval, see Table 1. The costate equation (50) reduces to
ψ T = ( 2 K ) 3 / 2 d w d K = 2 3 / 2 R ( K )
with R ( K ) = K 3 / 2 d w d K being a known and monotonely increasing function in K. Moreover, the function value R ( K ) is strictly positive for K > 0 . Since ψ T is constant, we have
K ( s ) R 1 2 3 / 2 ψ T = : K ptr ,
where R 1 stands for the inverse function of R. This means that, for the optimal solution, partial traction in the interior of the state domain is only possible at constant speed with K ( s ) K ptr . From equation (54) it follows that ψ T < 0 must hold.
Case 2: partial regenerative brake below speed limit. Likewise, if S i be an interval of partial regenerative brake in the interior of the state domain, then ψ ( s ) α and μ 3 ( s ) μ 4 ( s ) 0 will hold, and the costate equation now reads
ψ T = α ( 2 K ) 3 / 2 d w d K = 2 3 / 2 α R ( K ) .
Since ψ T is constant,
K ( s ) R 1 2 3 / 2 ψ T / α = : K pbr
holds, meaning that partial regenerative brake in the interior of the state domain is only possible at constant speed with K ( s ) K pbr . (In Khmelnitsky [1], the constants K ptr and K pbr are called K s and K b s , respectively.)
Case 3: partial traction on the speed limit. Let S i be an open interval of partial traction on the speed limit. Then ψ ( s ) 1 for s S i . Both u tr and μ 3 are unknown, but the equation of motion is entirely determined by the speed limit with v ( s ) = v max ( s ) and K ( s ) = K max ( s ) . Since v max is piecewise constant, and v can not be discontinuous, v, v max , K, and K max must be constant in the interval S i .
Case 4: partial regenerative brake on the speed limit. Likewise, if S i is an open interval of partial regenerative brake on the speed limit, then ψ ( s ) α for s S i . Both u br and μ 4 are unknown in this case, but the equation of motion is again entirely determined by the speed limit with v ( s ) v max and K ( s ) K max . Both v max and K max must be constant on S i , using the same argument as in Case 3.
Remark. In recent literature (not in Khmelnitsky), the constant speed driving regimes are often called cruising.

2.4. The algorithm of Khmelnitsky for fixed ψ T

In this section we are going to desribe the algorithm of finding a solution to the optimisation problem under the assumption that the Lagrange multiplier ψ T < 0 is a given number. Then, K ptr and K pbr are defined according to (54) and (56), respectively. Following Khmelnitsky [1], we define intervals with possible constant speed that correspond to the four cases studied in Section 2.3. In the case of constant speed, d K / d s = 0 holds, and the state differential equation (20) takes the form
u tr u br = d P d s + w .
Definition 1.
A sub-interval I of [ 0 , S ] with
0 d P ( s ) d s + w ( K ptr ) g tr ( K ptr ) and K ptr < K max ( s )
for all s I that cannot be extended, i.e., any enlargement of the interval would violate (58), is called a PT-interval. PT stands for partial traction.
Definition 2.
A sub-interval I of [ 0 , S ] with
g br ( K pbr ) d P ( s ) d s + w ( K pbr ) 0 and K pbr < K max ( s )
for all s I that cannot be extended is called a PB-interval. PB stands for partial regenerative brake.
Definition 3.
A sub-interval I of [ 0 , S ] with
0 d P ( s ) d s + w ( K max ( s ) ) g tr ( K max ( s ) ) and K max ( s ) = const .
for all s I that cannot be extended is called a PT-SL-interval. PT-SL stands for partial traction on speed limit.
Definition 4.
A sub-interval I of [ 0 , S ] with
g br ( K max ( s ) ) d P ( s ) d s + w ( K max ( s ) ) 0 and K max ( s ) = const .
for all s I that cannot be extended is called a PB-SL-interval. PB-SL stands for partial regenerative brake on speed limit.
Intervals of type PT, PB, PT-SL, and PB-SL are intervals where a constant speed motion of the optimum solution would be allowed by the control constraints. Those intervals are summarised under the name pcs-intervals, meaning ‘possible constant speed’. (In Khmelnitsky [1], PT is called minor grade, PB is called steep fall, PT-SL is called minor grade imitation, and PB-SL is called steep fall imitation.) Due to their definition, pcs-intervals will never intersect. In this paper, the start point at s = 0 , the pcs-intervals, and the stop point at s = S are summarised under the name ports. Additional speed limit ports might be introduced as will be explained later. All ports are numbered in the order of increasing s.
The optimal solution is found by connecting ports by trajectories of regular motion, namely full traction, coasting, or full regenerative brake. Below the speed limit, regular motion is governed by the differential equations
d K d s = d P d s + u w ,
d ψ d s = ψ d w d K + ψ T · ( 2 K ) 3 / 2 μ 3 d g tr d K μ 4 d g br d K .
with
u = g tr i f ψ > 1 0 i f α ψ 1 g br i f ψ < α
and μ 3 and μ 4 as defined in (51) and (), respectively. Equation (62) follows from (20) and Table 1, while () results from (50). Since u is discontinuous at ψ = α and ψ = 1 , the trajectory of K will have a kink point whenever ψ crosses these values.
It follows from the maximum principle that if the optimal solution touches the speed limit at a position s i , then h ( s i ) = 0 and η i might be positive by equation (45). This means for the first component of equation (43),
ψ ( s i + ) ψ ( s i ) = η i h E = η i h K = η i 0 ,
i.e., ψ might have a positive jump whenever the speed limit is attained. The height of this jump, however, is not defined in the maximum principle. It must be found by one-dimensional search.
If one tries to connect the ports A and B by a trajectory of regular motion, and this trajectory violates the speed limit, then a connection from A to B is not possible by this trajectory. In this case, one looks for a connection to the speed limit itself. If a trajectory is found that touches the speed limit in a single point, a new port is inserted there. It is called a speed limit touchpoint (SL-TP). Ports are now renumbered to be in the order of increasing s again. If the new speed limit touchpoint is located between the begin and the end position of a pcs-interval, then numbering is done such that the speed limit touchpoint comes first. If the speed limit touchpoint has number i, the ψ -value of the incoming trajectory is stored as ψ i . At the speed limit touchpoint, ψ is allowed to have a positive jump. A connection of ports in order to construct the optimal solution is only allowed in the direction of increasing port numbers.
When a trajectory of regular motion leaves a port, it is called a take-off, when it arrives at a port, it is called a landing. When trying to connect a port A with a port B by a trajectory of regular motion, there is always one of the variables s, K, and ψ not fixed at both take-off and landing. These variables can be adjusted to make the connection possible. The following cases of take-off and landing can exist:
  • (T1): take-off from the start point at s = 0 with K = 0 . The value of ψ is not fixed, but must be greater than 1 since full traction is applied.
  • (T2): take-off from interval of type PT with K = K ptr and ψ = 1 ± ϵ , with some small ϵ > 0 . Since d K / d s is discontinuous at ψ = 1 , the trajectories of K will leave in different directions depending on the choice of ψ slightly above or below 1, so both must be checked. The take-off position s is not fixed.
  • (T3): take-off from interval of type PB with K = K pbr and ψ = α ± ϵ , again with some small ϵ > 0 . The take-off position s is not fixed.
  • (T4): take-off from interior of a PT-SL interval with K = K max and ψ = 1 . The take-off position s is not fixed.
  • (T5): take-off from the end of a PT-SL interval with K = K max and ψ 1 , since a jump in ψ is allowed here.
  • (T6): take-off from interior of a PB-SL interval with K = K max and ψ = α . The take-off position s is not fixed.
  • (T7): take-off from the end of a PB-SL interval with K = K max and ψ α , since a jump in ψ is allowed here.
  • (T8): take-off from an SL-TP with number i: start with K = K max and ψ ψ i , since a jump in ψ is allowed here.
  • (L1): landing on PT interval with K = K ptr and ψ = 1 . The landing position s is not fixed.
  • (L2): landing on PB interval with K = K pbr and ψ = α . The landing position s is not fixed.
  • (L3): landing on the start of a PT-SL interval with K = K max and ψ < 1 . Then, a new speed limit touchpoint is inserted at the landing position, is connected with the PT-SL interval, and port renumbering is done such that the new speed limit touchpoint comes before the PT-SL interval. If the new speed limit touchpoint has number i, the ψ -value of the incoming trajectory is stored as ψ i .
  • (L4): landing on the start of a PT-SL interval with K = K max and ψ = 1 . Since any jump of ψ here would lead to full traction, it would violate the speed limit. Therefore, no new speed limit touchpoint is inserted.
  • (L5): landing on the interior of a PT-SL interval with K = K max and ψ = 1 . The landing position s is not fixed.
  • (L6): landing on the start of a PB-SL interval with K = K max and ψ < α . Then, a new speed limit touchpoint is inserted at the landing position, is connected with the PB-SL interval, and port renumbering is done such that the new speed limit touchpoint comes before the PB-SL interval. If the new speed limit touchpoint has number i, the ψ -value of the incoming trajectory is stored as ψ i .
  • (L7): landing on the start of a PB-SL interval with K = K max and ψ = α . Since any jump of ψ here would lead to coasting, it would violate the speed limit. Therefore, no new speed limit touchpoint is inserted.
  • (L8): landing on the interior of a PB-SL interval with K = K max and ψ = α . The landing position s is not fixed.
  • (L9): landing on an SL-TP with K = K max . The ψ -value of the incoming trajectory is stored, and ψ is allowed to jump here.
  • (L10): landing on the end point at s = S with K = 0 . The value of ψ is not fixed.
The construction of an optimal solution is best illustrated using a numerical example.
Example 1. A train is driven from Station A at s = 0 to Station B at s = S = 20 k m . The altitude is given by z = 40 m · sin ( s / km ) + 100 m , see Figure 2. Parameters are set according to Table 2.
In Example 1, no speed limit is assumed, and we consider the case ψ T = 1 . Equations (54) and (56) lead to K ptr = 536 m 2 / s 2 and K pbr = 665 m 2 / s 2 . Figure 3 shows the pcs-intervals calculated with equations (58) to (61), together with the kinetic energy of fastest motion. In our example, this leads to the exclusion of the first PT-intervall since it cannot be reached.
The first step in connecting ports would be the take-off from the start point at s = 0 . This is take-off case (T1). The trajectories of K and ψ are calculated according to equations (62) and (63). Different values of ψ at s = 0 lead to different trajectories that are shown in Figure 4 and Figure 5. When the ψ -trajectory crosses the value 1, full traction changes to coasting, and the corresponding K-trajectory has a kink point. When the ψ -trajectory crosses the value α = η tr η br = 0.7225 , coasting changes to full regenerative brake, and the corresponding K-trajectory has a kink point again. The blue line marks the only trajectory that would land on the PB-interval with number 2 (landing case (L2)). It is found by iterative bisection of the ψ -values at s = 0 .
Figure 6 and Figure 7 show the take-off from the PT-interval number 3. This is take-off case (T2). On the interval, the trajectories show an unstable behaviour with respect to ψ . If ψ is slightly above the value of 1, both K and ψ will move upwards. If ψ is slightly below the value of 1, both K and ψ will move downwards. Note that ψ takes off tangentially on the entire interval while K does so only from the ends of the interval. Close to the right end of the interval, the upwards moving K- and ψ -trajectories will soon turn downwards, and by that change the driving modus from full traction to coasting. This is a typical picture for take-off from PT- and PB-intervals.
Figure 8 and Figure 9 show the K- and ψ -trajectories when trying to connect intervals 4 and 8. The K-trajectories have a kink point when the driving modus changes from coasting to full traction or full regenerative brake, corresponding to ψ crossing the values 1 or α . Both the K- and ψ -trajectories are able to cross pcs-intervals, but the trajectory field often splits at pcs-intervals, as here at interval 7, where a shadowed region lies behind that cannot be reached by the trajectories.
Remark. There is great potential for speeding-up the algorithm by detection of those shadowed regions. For example, one can conclude immediately from Figure 4 that interval 3 cannot be reached from starting point 1, and from Figure 8 that interval 9 cannot be reached from interval 4. When solving the optimisation problem, it will always pay-off to invest into good statistics, showing which connections should be checked and which ones can safely be excluded.
Example 2. We consider Example 1, but now with a speed limit according to Table 3.
The multiplier ψ T is again set to 1 . Figure 10, Figure 11 and Figure 12 show that all types of pcs-intervals occur. The speed limit touchpoints 4, 10, and 14 have been inserted during the run of the algorithm.
Table 4 shows the possible port connections in Example 2.
It has been shown by Khmelnitsky [1] that there exists always exactly one connection from the start point (1) to the stop point (here 17). In Example 2, this is the connection 1 2 4 5 9 10 11 13 14 15 17 . The total time T ˜ required on this connection is not known a priori, but can be calculated from the K-curve. Since K = v 2 / 2 is the specific kinetic energy,
T ˜ = 0 S d s 2 K
holds.
Remark. Here, we distinguish between the scheduled time T, and the time T ˜ that is evaluated from the algorithm. The final goal of the algorithm is to match T ˜ to the prescribed T by variation of ψ T . This will be explained in Section 2.5.
Example 3. We consider Example 2, but now with ψ T = 0.3 . Trajectories of K and ψ are illustrated in Figure 13 and Figure 14, respectively.
The connections shown in Figure 13 and Figure 14 are listed in Table 5.

2.5. The algorithm of Khmelnitsky: variation of the multiplier ψ T

Khmelnitsky [1] has shown that the time T ˜ between two stops of the train is a strictly monotonously increasing function of the Lagrange multiplier ψ T . Moreover, ψ T < 0 and
lim ψ T 0 T ˜ =
hold. Therefore, whenever the scheduled time T is possible to be driven on a track section by a particular train, it can be approached by the optimisation algorithm by adjusting ψ T with iterative bisection. The complete algorithm for the minimum energy operation would read as follows:
Consider a track section with stops at s = 0 and s = S to be driven in a scheduled duration T.
  • Step 1: Calculate the fastest possible motion on the track section with purely regenerative braking. If the time T ˜ needed for that is larger than the scheduled time T, add mechanical braking and leave the algorithm. If T ˜ < T , choose an arbitrary negative value for ψ T and proceed the algorithm with Step 2.
  • Step 2: Calculate K ptr by (54), K pbr by (56), and evaluate the pcs-intervals PT, PB, PT-SL, PB-SL by (58)–(61). Skip pcs-intervals that cannot be reached by fastest motion. Numerate the ports, i.e., the remaining pcs-intervals, the start at s = 0 , and the stop at s = S , in the order of ascending s.
  • Step 3: Try to connect ports by regular motion with equations (62) and (). Use one of the take-off cases (T1)–(T8) and one of the landing cases (L1)–(L10). Add speed limit touchpoints if necessary, according to the instructions given above. Step 3 is complete when a connection from the start point at s = 0 to the stop point s = S has been found.
  • Step 4: Calculate T ˜ according to (66). If T ˜ is sufficiently close to T, the algorithm is successfully completed. If not, adjust ψ T and proceed with Step 2.
If n is the number of ports then a maximum of n 2 = ( n 2 n ) / 2 possible connections has to be checked, unless the start-stop connection is found earlier. This means that the number of possible connections grows quadratically with n. The algorithm can be seen as a search tree, with port connections being the branches of the tree that need to be checked. Therefore, it is crucial to follow a clever search strategy to keep computing time at an acceptable level. We mention four important measures that dramatically shortened computing time when the code was developed:
Parallel path exclusion. Khmelnitsky [1] has shown that any two ports can only be connected by a at most one path. Therefore, it is wise to exclude all parallel paths in the search tree. For example, if a connection 1 2 3 has been established, the parallel direct link 1 3 is not possible and need not be checked.
Start from the treetop. When choosing the next connection to check for in Step 3, start at the highest port number that is connected to port 1 and not a known dead end. This strategy usually leads to an early discovery of the start-stop connection.
Look out for shadowed ports. If ports lie in shadow, exclude impossible connections. See the remark in the discussion of Example 1 in Section 2.3.
Estimate ψ T by interpolation. The adjustment of ψ T in Step 4 can be sped up using interpolation techniques.
Example 4. Example 4 is equal to Examples 2 and 3, except that ψ T is not prescribed, and T is set to 16 minutes.
Applying the algorithm, ψ T converges to 0.5237 . The final K- and ψ -trajectories are shown in Figure 15 and Figure 16.
Connections and type of take-off and landing are given for the final value ψ T = 0.5237 in Table 6.
In Figure 17 and Figure 18, speed and electric energy are shown for Example 4 at ψ T = 0.5237 .

3. Results and discussion

3.1. Estimating the potential of efficient train driving

In order to compare the optimal driving performance to the energy demand of trains in real-life operation, a study has been carried out covering a total of 100 regional train runs on four track sections in Germany. The trains were drawn by electric locomotives equipped with regenerative braking. Time, velocity, energy supply from the catenary, and energy return by regenerative braking have been recorded. Trains with five, six, or seven coaches were included in the study. Additional energy requirements E add , as for air-conditioning, ventilation, air compression, etc. have been estimated to an average value of 254 kW per train from the measured data. Train runs affected by speed restriction due to signalling in the interior of the track section have been excluded.
In the study it would not be approapriate to compare the energy demand of a particular train run to the energy optimum based on the duration given by the timetable. Trains that are delayed need to drive faster, so the shorter time they have available needs to be considered in the optimum calculation. On the other hand, trains that arrive early are seen to unnecessarily waste energy by driving too fast. Therefore, for each train run a reference time is considered that allowes a fair comparison to the optimum. Let
  • t start TT be the departure time according to timetable,
  • t stop TT be the arrival time according to timetable,
  • t start rec be the recorded departure time, and
  • t stop rec be the recorded arrival time of the train run.
Then, the reference time duration is defined by T = max ( t stop TT , t stop rec ) t start rec , and the energy optimum is calculated with respect to this reference time.
A total energy demand is defined by E tot = E tr E br + E add , where E tr is electric energy used for traction, E br is electric energy returned by regenerative braking, and E add = 254 k W is the above mentioned additional energy requirement. Figure 19, Figure 20, Figure 21 and Figure 22 display the total energy E tot over the reference time T. The stars in the figures are recorded measurements, while the lines indicate the calculated energy optimum.
Result.Figure 19, Figure 20, Figure 21 and Figure 22 show that the recorded total energy spreads over a rather wide range. A factor of approximately two lies between minimum and maximum energy consumption on each track section, even for cases with quite similar reference time. Contrary to the expectation, the recorded values do not show an increase of energy demand with train length, but this might be due to the small size of the sample and the wide spread of values. For many train runs, the measured energy is far above the corresponding optimum line, in some extremal cases by a factor of three.
Remark. In one case shown in Figure 19, a seven-coach train was in fact better than the optimum. This is not a contradiction to the optimality property. As it has been pointed out, optimal control theory ensures that the algorithm returnes the energy minimum, but this only holds for the given train motion model. This model is only a simplification, and some assumptions might not have been met during the train drive experiment. For example, wind speed is neither measured nor considered in the model, and the additional energy requirement E add is only averaged for a lack of individual measurements.
If for the i-th train run the recorded total energy is denoted by E tot , i rec , and the corresponding optimum energy is E tot , i opt , then the energy ratio η i = E tot , i opt / E tot , i rec can be evaluated. Within the concept of the Physical Optimum (PhO) [45,46], η i is called a PhO factor that quantitatively evaluates the efficiency with respect to a feasible most efficient reference process. This energy ratio is shown in Figure 23 for each of the 100 train runs of the study. A total energy ratio η tot = i E tot , i opt / i E tot , i rec = 0.627 is obtained by summing over all train runs studied. This can be seen as an estimate of the potential that energy optimised driving would have. Based on the results of the present study, an amount of 37.3 percent of energy could be saved by following an energy minimising driving strategy. One should, however, bear in mind two limitations: First, due to model simplifications, the computed driving strategy might differ from the real optimum under the current conditions. And second, even if an optimum strategy is displayed to the driver by an assistance system, it can only be approached, and interactions with other traffic will sometimes not allow to exactly follow the suggestions. However, since the theoretical optimum turned out to be substantially lower than the measured energy demand in operation, we think that there is still a large potential for energy saving by optimal control based assistance for train drivers.

4. Conclusion

Within this paper, energy-efficient train driving has been studied for the case of an exclusive usage of the regenerative brake in electric trains. The code opTop written at Fraunhofer Institute Magdeburg is based on the algorithm of Khmelnitsky that constructs the unique minimum energy solution. A derivation of the statements given by Khmelnitsky from a more general formulation of Pontryagin’s maximum principle is presented. In addition to the theory in Khmelnitsky’s article, a complete list of switching cases has been provided and illustrated by a number of numerical examples. A comparison to energy consumption data in real operation showed that the energy minimising strategy was able to save, on average, about 37 % of energy. Extensions of the code to include mechanical braking, non-zero speed boundary conditions, and dynamic response to train delays are subject to further research.

Author Contributions

Writing—original draft preparation, W.H.; software, W.H.; project administration, M.R. and T.B.; funding acquisition, M.R. and T.B. All authors have read and agreed to the published version of the manuscript.

Funding

The authors gratefully acknowledge external funding by the

Data Availability Statement

External data that have been processed cannot be provided nor referenced for confidentiality reasons.

Acknowledgments

The authors would like to thank the team of traveltainer GmbH & Co. KG, Altenbeken, Germany for their valuable collaboration and support during this research project.

Conflicts of Interest

The authors declare no conflict of interest.

Abbreviations

The following abbreviations are used in this manuscript:
DDP discrete dynamic programming
EETC energy-efficient train control
KKT Karush-Kuhn-Tucker (conditions)
NLP nonlinear programming
PB partial regenerative braking
PB-SL partial regenerative braking at speed limit
pcs-interval possible constant speed interval
PMP Pontryagin maximum principle
PT partial traction
PT-SL partial traction at speed limit
SL-TP speed limit touchpoint
SNCF Société nationale des chemins de fer français (French national railways)
TGV train à grande vitesse (French high speed train)

References

  1. Khmelnitsky, E. On an Optimal Control Problem of Train Operation. IEEE Transaction on Automatic Control 2000, 45, 1257–1266. [Google Scholar] [CrossRef]
  2. Strobel, H.; Horn, P.; Kosemund, M. A contribution to optimum computer-aided control of train operation. In Proceedings of the Proceedings of the 2nd FAC/IFIP/IFORS symposium on traffic control and transportation systems, 1974, pp.; pp. 377–387.
  3. Howlett, P.G.; Pudney, P.J. Energy-efficient train control; Springer: London, 1995. [Google Scholar]
  4. Franke, R.; Terwiesch, P.; Meyer, M. An algorithm for the optimal control of the driving of trains. In Proceedings of the Proceedings of the 39th IEEE conference on decision and control, 2000, Vol.; p. 3.
  5. Albrecht, T. Ein Beitrag zur Nutzbarmachung genetischer Algorithmen für die optimale Steuerung und Planung eines flexiblen Stadtschnellbahnbetriebes (A contribution to the utilization of Genetic Algorithms for optimal control and planning of a flexible urban rail transport). PhD thesis, 2005. [Google Scholar]
  6. Howlett, P.G.; Pudney, P.J.; Vu, X. Freightmiser: an energy-efficient application of the train control problem. In Proceedings of the Proceedings of the 30th conference of australian institutes of transport research (caitr); 2008. [Google Scholar]
  7. Aradi, S.; Becsi, T.; Gaspar, P. A predictive optimization method for energy-optimal speed profile generation for trains. In Proceedings of the Proceedings of the 2013 IEEE 14th international symposium on computational intelligence and informatics (CINTI); 2013. [Google Scholar]
  8. DAS on-board algorithms. Technical Report D6.3, ON-TIME: Optimal Networks for Train Integration Management across Europe, 2014.
  9. Albrecht, A.R.; Howlett, P.G.; Pudney, P.J.; Vu, X.; Zhou, P. Energy-efficient train control: The two-train separation problem on level track. Journal of Rail Transport Planning & Management 2015, 5, 163–182. [Google Scholar]
  10. Scheepmaker, G.; Goverde, R.; Kroon, L. Review of energy-efficient train control and timetabling. European Journal of Operational Research 2017, 257, 355–376. [Google Scholar] [CrossRef]
  11. Pontryagin, L.; Boltyanskii, V.; Gamkrelidze, R.; Mishchenko, E. Matematicheskaya teoriya optimal’nykh prozessov; Fizmatgiz: Moscow, 1961. [Google Scholar]
  12. Pontryagin, L.; Boltyanskii, V.; Gamkrelidze, R.; Mishchenko, E. The Mathematical Theory of Optimal Processes; John Wiley and Sons (Interscience Publishers): New York, 1962. [Google Scholar]
  13. Ichikawa, K. Application of optimization theory for bounded state variable problems to the operation of train. Bulletin of the Japan Society of Mechanical Engineers 1968, 11, 857–865. [Google Scholar] [CrossRef]
  14. Milroy, I.P. Aspects of automatic train control. PhD thesis, Loughborough University, Loughborough, UK, 1980. PhD thesis, Loughborough University, Loughborough, UK, 1980. [Google Scholar]
  15. Howlett, P.G.; Milroy, I.P.; Pudney, P.J. Energy-efficient train control. Control Engineering Practice 1994, 2, 193–200. [Google Scholar] [CrossRef]
  16. Liu, R.R.; Golovitcher, I.M. Energy-efficient operation of rail vehicles. Transportation Research Part A: Policy and Practice 2003, 37, 917–932. [Google Scholar] [CrossRef]
  17. Albrecht, A.R.; Howlett, P.G.; Pudney, P.J.; Vu. Energy-efficient train control: From local convexity to global optimization and uniqueness. Automatica 2003, 49, 3072–3078. [Google Scholar] [CrossRef]
  18. Asnis, I.A.; Dmitruk, A.V.; Osmolovskii, N.P. Solution of the problem of the energetically optimal control of the motion of a train by the maximum principle. USSR Computational Mathematics and Mathematical Physics 1985, 25, 37–44. [Google Scholar] [CrossRef]
  19. Ying, P.; Zeng, X.; Song, H.; Shen, T.; Yuan, T. Energy-efficient train operation with steep track and speed limits: A novel Pontryagin’s maximum principle-based approach for adjoint variable discontinuity cases. IET Intelligent Transport Systems 2021, 15, 1183–1202. [Google Scholar] [CrossRef]
  20. Ying, P.; Zeng, X.; Shen, T.; Wang, Y.; Ma, Z.; Wu, Y. Partial Train Speed Trajectory Optimization. In Proceedings of the 2022 IEEE 25th International Conference on Intelligent Transportation Systems (ITSC), Macau, China; 2022. [Google Scholar]
  21. Baranov, L.A.; Meleshin, I.S.; Chin’. Optimal control of a subway train with regard to the criteria of minimum energy consumption. Russian Electrical Engineering 2011, 82, 405–410. [Google Scholar] [CrossRef]
  22. Lu, S.; Weston, P.; Hillmansen, S.; Gooi, H.B. andRoberts, C. Increasing the regenerative braking energy for railway vehicles. IEEE Transactions on Intelligent Transportation Systems 2014, 15, 2506–2515. [Google Scholar] [CrossRef]
  23. Zhou, Y.; Bai, Y.; Li, J.; Mao, B.; Li, T. Integrated optimization on train control and timetable to minimize net energy consumption of metro lines. Journal of Advanced Transportation 2018, 1–19. [Google Scholar] [CrossRef]
  24. Fernández-Rodríguez, A.; Fernández-Cardador, A.; Cucala, A. Real time eco-driving of high speed trains by simulation-based dynamic multi-objective optimization. Simulation Modelling Practice and Theory 2018, 84, 50–68. [Google Scholar] [CrossRef]
  25. Scheepmaker, G.; Goverde, R. Energy-efficient train control using nonlinear bounded regenerative braking. Transportation Research Part C: Emerging Technologies 2020, 121, 102852. [Google Scholar] [CrossRef]
  26. Rao, A.; Benson, D.; Darby, C.; Patterson, M.; Francolin, C.; Sanders, I.; Huntington, G. Algorithm 902: GPOPS, a MATLAB software for solving multiple-phase optimal control problems using the Gauss pseudospectral method. ACM Transactions an Mathematical Software 2010, 37, 22:1–22:39. [Google Scholar] [CrossRef]
  27. Ghaviha, N.; Bohlin, M.; Holmberg, C.; Dahlquist, E.; Skoglund, R.; Jonasson, D. A driver advisory system with dynamic losses for passenger electric multiple units. Transportation Research Part C: Emerging Technologies 2017, 85, 111–130. [Google Scholar] [CrossRef]
  28. Kouzoupis, D.; Pendharkar, I.; Frey, J.; Diehl, M.; Corman, F. Direct multiple shooting for computationally efficient train trajectory optimization. Transportation Research Part C: Emerging Technologies 2023, 152, 104170. [Google Scholar] [CrossRef]
  29. Feng, M.; Huang, Y.; Lu, S. Eco-driving Strategy Optimization for High-speed Railways Considering Dynamic Traction System Efficiency. IEEE Tranactions on Transportation Electrification 2023, 1–11. [Google Scholar] [CrossRef]
  30. Szkopiński, J.; Kochan, A. Energy Efficiency and Smooth Running of a Train on the Route While Approaching Another Train. Energies 2021, 14, 7593. [Google Scholar] [CrossRef]
  31. Wang, Y.; Schutter, B.; van den Boom, T.; Ning, B. Optimal trajectory planning for trains - A pseudospectral method and a mixed integer linear programming approach. Transportation Research Part C: Emerging Technologies 2013, 29, 97–114. [Google Scholar] [CrossRef]
  32. Ross, I.; Karpenko, M. A review of pseudospectral optimal control: from theory to flight. Annual Reviews in Control 2012, 36, 182–197. [Google Scholar] [CrossRef]
  33. Ross, I. A Primer on Pontryagin’s Principle in Optimal Control; Collegiate Publishers: San Fransisco, USA, 2015. [Google Scholar]
  34. Becerra, V. Solving complex optimal control problems at no cost with psopt. In Proceedings of the 2010 IEEE International Symposium on Computer-Aided Control System Design, Yokohama, Japan; 2010; pp. 1391–1396. [Google Scholar]
  35. Hartl, R.; Sethi, S.; Vickson, R. A survey of the maximum principles for optimal control problems with state constraints. SIAM Review 1995, 37, 181–218. [Google Scholar] [CrossRef]
  36. Ihme, J. Schienenfahrzeugtechnik; Springer Fachmedien: Wiesbaden, 2016. [Google Scholar]
  37. Fassbinder, S. Wie Energie-effizient ist der Bahnverkehr wirklich? Technical report, Deutsches Kupferinstitut, Düsseldorf, 2010. Available online: https://de.slideshare.net/LeonardoENERGYDeutschland/wie-energieeffizient-ist-der- bahnverkehr-wirklich.
  38. Curtius, E.; Kniffler, A. Neue Erkenntnisse über Haftung zwischen Treibrad und Schiene. Elektrische Bahnen 1944, 1950, 1/2, 3/4, 9. /2.
  39. Schlecht, B. Maschinenelemente 2; Pearson Studium: Munich, 2010. [Google Scholar]
  40. Wie grün“ ist der Schienenverkehr? Technical report, Deutsches Kupferinstitut, Düsseldorf. 2020. Available online: https://www.kupferinstitut.de/wp-content/uploads/2020/02/BahnEffizienz.pdf.
  41. Frilli, A.; Meli, E.; Nocciolini, D.; Pugi, L.; Rindi, A. Energetic optimization of regenerative braking for high speed railway systems. Energy Conversion and Management 2016, 129, 200–215. [Google Scholar] [CrossRef]
  42. Voß, G.; Gackenholz, L.; Wiebels, R. Eine neue Formel (Hannoversche Formel) zur Bestimmung des Luftwiderstandes spurgebundener Fahrzeuge. Glasers Annalen – Zeitschrift für Eisenbahnwesen und Verkehrstechnik 1972, 96, 166–171. [Google Scholar]
  43. Karush, W. Minima of Functions of Several Variables with Inequalities as Side Constraints. Master’s thesis, University of Chicago, Chicago, Illinois, 1939. [Google Scholar]
  44. Kuhn, H.; Tucker, A. Nonlinear Programming. In Proceedings of the Proceedings of 2nd Berkeley Symposium, Berkeley, 1951; pp. 481–492.
  45. Volta, D. Das Physikalische Optimum als Basis von Systematiken zur Steigerung der Energie- und Stoffeffizienz von Produktionsprozessen. PhD thesis, Clausthal University of Technology, 2014. [Google Scholar]
  46. Eggers, N.; Böttger, J.; Kerpen, L.; Sankol, B.; Birth, T. Refining VDI guideline 4663 to evaluate the efficiency of a power-to-gas process by employing limit-oriented indicators. Energy Efficiency 2021, 14, 1–15. [Google Scholar] [CrossRef]
Figure 1. Maximum force for traction and regenerative brake.
Figure 1. Maximum force for traction and regenerative brake.
Preprints 82118 g001
Figure 2. Track altitude in Example 1.
Figure 2. Track altitude in Example 1.
Preprints 82118 g002
Figure 3. Example 1: pcs-intervals PT and PB, fastest motion, and numbering.
Figure 3. Example 1: pcs-intervals PT and PB, fastest motion, and numbering.
Preprints 82118 g003
Figure 4. Example 1: K-trajectories starting from s = 0 (port 1) for various values of ψ at s = 0 . Trajectories change from full traction to coasting and then to full regenerative braking. The blue trajectory connects port 1 with port 2.
Figure 4. Example 1: K-trajectories starting from s = 0 (port 1) for various values of ψ at s = 0 . Trajectories change from full traction to coasting and then to full regenerative braking. The blue trajectory connects port 1 with port 2.
Preprints 82118 g004
Figure 5. Example 1: ψ -trajectories starting from s = 0 (port 1) for various values of ψ at s = 0 . Trajectories change from full traction to coasting and then to full regenerative braking. The blue trajectory connects port 1 with port 2.
Figure 5. Example 1: ψ -trajectories starting from s = 0 (port 1) for various values of ψ at s = 0 . Trajectories change from full traction to coasting and then to full regenerative braking. The blue trajectory connects port 1 with port 2.
Preprints 82118 g005
Figure 6. Example 1: K-trajectories starting from PT-interval with number 3.
Figure 6. Example 1: K-trajectories starting from PT-interval with number 3.
Preprints 82118 g006
Figure 7. Example 1: ψ -trajectories starting from PT-interval with number 3.
Figure 7. Example 1: ψ -trajectories starting from PT-interval with number 3.
Preprints 82118 g007
Figure 8. Example 1: K-trajectories starting from interval 4, heading for interval 8.
Figure 8. Example 1: K-trajectories starting from interval 4, heading for interval 8.
Preprints 82118 g008
Figure 9. Example 1: ψ -trajectories starting from interval 4, heading for interval 8.
Figure 9. Example 1: ψ -trajectories starting from interval 4, heading for interval 8.
Preprints 82118 g009
Figure 10. Example 2: K-trajectories.
Figure 10. Example 2: K-trajectories.
Preprints 82118 g010
Figure 11. Example 2: ψ -trajectories.
Figure 11. Example 2: ψ -trajectories.
Preprints 82118 g011
Figure 12. Example 2: zoom into ψ -trajectories.
Figure 12. Example 2: zoom into ψ -trajectories.
Preprints 82118 g012
Figure 13. Example 3: K-trajectories.
Figure 13. Example 3: K-trajectories.
Preprints 82118 g013
Figure 14. Example 3: ψ -trajectories.
Figure 14. Example 3: ψ -trajectories.
Preprints 82118 g014
Figure 15. Example 4: K-trajectories.
Figure 15. Example 4: K-trajectories.
Preprints 82118 g015
Figure 16. Example 4: ψ -trajectories.
Figure 16. Example 4: ψ -trajectories.
Preprints 82118 g016
Figure 17. Example 4: train speed.
Figure 17. Example 4: train speed.
Preprints 82118 g017
Figure 18. Example 4: net electric energy.
Figure 18. Example 4: net electric energy.
Preprints 82118 g018
Figure 19. Total energy demand E tot over reference time T for track section 1.
Figure 19. Total energy demand E tot over reference time T for track section 1.
Preprints 82118 g019
Figure 20. Total energy demand E tot over reference time T for track section 2.
Figure 20. Total energy demand E tot over reference time T for track section 2.
Preprints 82118 g020
Figure 21. Total energy demand E tot over reference time T for track section 3.
Figure 21. Total energy demand E tot over reference time T for track section 3.
Preprints 82118 g021
Figure 22. Total energy demand E tot over reference time T for track section 4.
Figure 22. Total energy demand E tot over reference time T for track section 4.
Preprints 82118 g022
Figure 23. Energy ratio η i for each run, shown in descending order.
Figure 23. Energy ratio η i for each run, shown in descending order.
Preprints 82118 g023
Table 1. Case distinction for ψ . ‘cont.’ means that the continuity of the Lagrange multipliers has been exploited. The corresponding value is only valid in the interior of the state domain, i.e., when v < v max .
Table 1. Case distinction for ψ . ‘cont.’ means that the continuity of the Lagrange multipliers has been exploited. The corresponding value is only valid in the interior of the state domain, i.e., when v < v max .
case u tr = u br = μ 1 = μ 2 = μ 3 = μ 4 =
ψ > 1 g tr 0 0 ψ α ψ 1 0
(full traction) (38) (38) (41) (49) (48) (41)
ψ = 1 0 0 1 α 0 0
(partial traction) (38) cont. (49) cont. (41)
α < ψ < 1 0 0 1 ψ ψ α 0 0
(coasting) (38) (38) (48) (41) (41) (41)
ψ = α 0 1 α 0 0 0
(partial reg. brake) (38) (48) cont. (41) cont.
ψ < α 0 g br 1 ψ 0 0 α ψ
(full reg. brake) (38) (38) (48) (41) (41) (49)
Table 2. Parameters of Example 1.
Table 2. Parameters of Example 1.
parameter symbol value
number of coaches n 6
mass of locomotive m loc 84 t
total mass of train m 414 t
engine efficiency η tr 0.85
efficiency of regenerative brake η br 0.85
max. mechanical power for traction P tr , mech 5.6 M W
max. mechanical power for regenerative brake P br , mech 5.6 M W
max. force of regenerative brake F br , lim 240 k N
air drag coefficient for locomotive c w , loc 0.26
air drag coefficient for first coach c w , first 0.13
air drag coefficient for middle coaches c w , middle 0.10
air drag coefficient for last coach c w , last 0.23
rolling friction coefficient c roll 0.0015
coefficient accounting for rotating masses c rot 1.08
Table 3. Speed limits in Example 2.
Table 3. Speed limits in Example 2.
from position s [km] to position s [km] max. speed v max [km/h]
0 5.5 160
5.5 7.0 110
7.0 9.6 150
9.6 12.0 105
12.0 20.0 140
Table 4. Port connections in Example 2.
Table 4. Port connections in Example 2.
port connection take-off and landing type remark
1 2 (T1), (L2) K has kink point when full traction changes to coasting.
2 4 5 (T3), (L3) New SL-TP 4 included, connected to 5, ψ jumps at 4.
5 6 (T5), (L5) ψ jumps at end of 5.
5 9 (T5), (L4) ψ jumps at end of 5. Since ψ = 1 at the beginning of interval 9, no new speed limit touchpoint is included.
6 7 (T5), (L1) ψ jumps at end of 6.
9 10 11 (T7), (L3) A connection with constant speed, but with a jump of ψ at new SL-TP 10, which is connected to 11.
11 12 (T5), (L1) ψ jumps at end of 11.
11 13 (T5), (L2) ψ jumps at end of 11. Same K-trajectory as 11 12 at the beginning, but then change to coasting.
13 14 (T3), (L9) Landing at newly included SL-TP 14.
14 15 (T8), (L1) Take-off from SL-TP 14 with ψ -jump.
15 16 (T2), (L1) Full traction.
15 17 (T2), (L10) Coasting to stop.
Table 5. Port connections in Example 3.
Table 5. Port connections in Example 3.
port connection take-off and landing type remark
1 2 (T1), (L1)
1 3 (T1), (L2) K-trajectories of 1 2 and 1 3 coincide at the beginning.
3 4 (T3), (L9) New SL-TP 4 is inside the s-range of interval 5. SL-TP 4 needs to come first in numbering.
4 5 (T8), (L1) ψ jumps at SL-TP 4.
5 6 (T2), (L2)
6 7 (T3), (L9) SL-TP 7 is inserted, lying exactly between intervals 6 and 8.
7 8 (T8), (L1) ψ jumps at SL-TP 7.
8 9 (T2), (L2)
8 10 (T2), (L1) Trajectory is close to speed limit but does not touch.
10 11 (T2), (L10) Coasting to stop.
Table 6. Port connections in Example 4 at final value ψ T = 0.5237 .
Table 6. Port connections in Example 4 at final value ψ T = 0.5237 .
port connection take-off and landing type remark
1 2 (T1), (L2)
2 3 (T3), (L9) New SL-TP 3 is inside the s-range of interval 4. The SL-TP 3 needs to come first in numbering.
3 4 (T8), (L1) ψ jumps at SL-TP 3.
4 5 (T2), (L2)
4 6 (T2), (L7) Landing on start of interval 6 with ψ = α . Any jump of ψ here would lead to coasting and violate speed limit. Therefore, no new speed limit touchpoint is inserted.
6 7 (T7), (L1) ψ jumps at end of PB-SL 6.
7 8 (T2), (L2)
7 9 (T2), (L9) New SL-TP 9 is exactly between intervals 8 and 10.
9 10 (T8), (L1) Take-off from SL-TP 9.
10 11 (T2), (L10) Coasting and finally full regenerative braking to stop.
Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.
Copyright: This open access article is published under a Creative Commons CC BY 4.0 license, which permit the free download, distribution, and reuse, provided that the author and preprint are cited in any reuse.
Prerpints.org logo

Preprints.org is a free preprint server supported by MDPI in Basel, Switzerland.

Subscribe

© 2024 MDPI (Basel, Switzerland) unless otherwise stated