1 Introduction

Driven by the objectives of lower emissions and lower operational costs, railway train energy simulation is a popular research topic [1,2,3]. Basically, there are two types of train models that are used to calculate train energy consumptions. The first type considers the train as a whole, in which inter-vehicle movements are neglected [2,3,4,5,6,7,8,9,10]. The second type uses models that are commonly used for longitudinal train dynamics (LTD) [11, 12]; specifically, each vehicle of the train is modelled as a rigid body with a single degree of freedom (DOF) [13,14,15,16,17,18,19]. Both types of models use the same or similar procedures to calculate train energy consumptions. These procedures determine locomotive traction torques or forces first, and then work out how much energy is consumed by the locomotives.

Traction force is the tangential (mostly in the longitudinal direction) friction force in wheel–rail contact; it is directly related to the coefficient of adhesion which is defined as the ratio of realised tangential force over the normal force. During train operations, traction forces are influenced by a number of factors, for example (1) vehicle dynamics that can vary wheel–rail vertical force, (2) wheel–rail adhesion characteristics that directly influence the resulted traction forces, and (3) traction control that changes traction torque to avoid excessive wheel slip and/or slide. However, none of these factors is able to be considered using the two conventional types of models. To better study train energy issues, this paper develops advanced train models that consider locomotive adhesion and LTD.

Although the modelling and simulation methods presented in this paper can also be adapted to be used for passenger trains, this paper focuses on train energy issues in freight trains including heavy haul trains. Compared to short passenger trains, the energy issues in freight trains have two differences. First, freight trains are usually significantly heavier than passenger trains, being 10–100 times greater. Considering these significant mass differences, the energy possessed by freight trains can also be significantly larger than by passenger trains. Second, many freight trains, especially heavy haul trains, run between mines and ports. There is an evident terrain feature between mines and ports. Mines are usually located in mountainous areas with relatively higher altitudes, whilst ports are all located seaside with almost zero altitude. Considering that most heavy haul trains return to mines empty, theoretically, there is a large amount of energy that can be recovered during the loaded trip via technologies such as dynamic brake (DB) and energy storage systems (ESSs) [3, 7, 9, 11, 12]. Studies that focus on heavy haul trains therefore have the potential to be of great significance.

The rest of this paper is organised as follows. Section 2 introduces the basics of LTD that will be used for subsequent modelling works. Section 3 presents different locomotive models that are used in this paper. Section 4 presents the simulation strategies, whilst Sect. 5 presents simulation case studies and resulting discussion.

2 LTD basics

LTD is a branch of vehicle system dynamics which simplifies each vehicle of a train as a single DOF rigid body. Only the longitudinal DOF of a vehicle is considered. LTD usually studies inter-vehicle movements of individual vehicles as well as the movements of the train as a whole. Force components that are considered in the system usually include in-train force (also called coupler force), traction force, air brake force, dynamic brake force, propulsion resistance force, grade force, and curving resistance force. The basic equations of motion can be expressed as

$$ ma = F_{{{\text{c}}1}} + F_{{{\text{c}}2}} + F_{\text{t}} + F_{\text{db}} + F_{\text{ab}} + F_{\text{pr}} + F_{\text{cr}} + F_{\text{G}} , $$
(1)

where \( m \) and \( a \) are the mass and acceleration of the vehicle, respectively; \( F_{{{\text{c}}1}} \) is the front coupler force; \( F_{{{\text{c}}2}} \) is the rear coupler force; \( F_{\text{t}} \) is the traction force; \( F_{\text{db}} \) is the DB force; \( F_{\text{ab}} \) is the air brake force; \( F_{\text{pr}} \) is the propulsion resistance force; \( F_{\text{cr}} \) is the curving resistance force; and \( F_{\text{G}} \) is the grade force. LTD is a popular research topic which has been discussed in many publications; more detailed information of LTD can be found in References [18, 19].

LTD is used in this paper for energy consumption simulations. Due to the simplifications of simulated wagons, i.e. only one DOF is considered for each wagon, the accuracy of energy calculations will also be influenced. For example, wheel–rail friction that has implications for train energy consumptions cannot be accurately simulated in LTD. However, the single DOF simplification in LTD is widely acknowledged and well accepted for train energy simulations in academia and industry [11,12,13,14,15,16,17].

3 Locomotive models

To better investigate the implications of locomotive wheel–rail adhesions for train energy calculations, three locomotive models are used in this paper: (1) a model commonly used for LTD simulations, which can be regarded as a one-dimensional model; (2) a three-dimensional (3D) model developed using a commercial software called GENSYS [20]; and (3) a two-dimensional (2D) in-house model developed in FORTRAN language.

3.1 1D locomotive model

A locomotive model that is commonly used for conventional LTD simulations is shown in Fig. 1. The model has considered (1) notch variable traction performance, (2) speed variable traction performance, (3) constant adhesion limit, and (4) gradual traction motor response. During simulations, traction notches and locomotive speeds are used to interpolate traction forces from locomotive traction performance curves. These forces are then modified by traction motor response limits and track conditions. Basically, this 1D locomotive model is based on look-up tables with empirical limits. The adhesion limit is an empirical constant parameter rather than parameters that are delivered from dynamics simulations. Wheel–rail contact is not considered at all. More information of this type of locomotive model can be found in Ref. [21].

Fig. 1
figure 1

1D locomotive model used for conventional LTD simulations [21]

Full size image

3.2 3D locomotive models

A 3D locomotive model was developed using the commercial software called GENSYS. The simulated locomotive is a diesel-electric heavy haul locomotive with a C0-C0 wheelset arrangement; it has an axle load of 22.3 t and bogie spacing of 14.2 m. Detailed model parameters are listed in “Appendix”. The locomotive model considers six DOFs for all major components: one car-body, two bogie frames, and six wheelsets. The model considers springs, dampers and bumpstops for primary suspensions, and the same plus bogie pivot pins for secondary suspensions.

The wheel–rail contact model has profiles shown in Fig. 2; the geometrical model allows a maximum of three points of contact. Normal forces in the contact model were calculated using a spring–damper element, whilst the tangential forces were calculated using the modified FASTSIM algorithm. The coefficient of friction in the modified FASTSIM is expressed as [22]

$$ \mu = \mu_{\text{s}} \left( {\left( {1 - A} \right){\text{e}}^{ - B\omega } + A} \right) , $$
(2)

where \( \mu \) is the actual coefficient of friction; \( \mu_{\text{s}} \) is the static coefficient of friction; \( A \) is ratio of the static coefficient of friction to the kinetic coefficient of friction; \( B \) is the coefficient of exponential friction decrease; and \( \omega \) is the creepage or slip. In this model, Kalker’s original constant reduction factor is modified to be a variable reduction factor that is expressed as [22]

$$ k = k_{0} \left( {\alpha_{ \inf } + \frac{{1 - \alpha_{ \inf } }}{1 + \beta \varepsilon }} \right) , $$
(3)

where \( k \) is the variable reduction factor; \( k_{0} \) is Kalker’s original reduction factor; \( \alpha_{ \inf } \) is the fraction of the initial value of Kalker’s reduction factor at creep values approaching infinity; \( \beta \) is a non-dimensional parameter related to the decrease in the slip area; and \( \varepsilon \) is a parameter describing the gradient of the tangential stress in the stress distribution transferred to hemisphere.

Fig. 2
figure 2

Wheel–rail profiles

Full size image

In the locomotive model, the traction process of this model is controlled by using a proportional–integral adhesion controller [23]; the slip limit for adhesion was set to 7%. One example of adhesion characteristics under the traction control is shown in Fig. 3.

Fig. 3
figure 3

Adhesion characteristics

Full size image

3.3 2D locomotive models

An in-house 2D locomotive model was developed in FORTRAN programming language. The 2D locomotive model is based on the same locomotive for the 3D model, and is basically a simplified version of the latter. The 2D model has 27 DOFs considering three DOFs (longitudinal translation, bounce, and pitch) for each of the bodies considered in the model: six wheelsets, two bogie frames, and one car-body. The wheel–rail contact model of the 2D locomotive model uses the same normal force and tangential force models that were used for the 3D locomotive model but assumes that the wheelsets are always placed in the centre of the track. In this way, the contact geometry of the wheel–rail interface will be constant throughout the simulations. One of the limitations of the simplification is that this model cannot consider curved track. When used for curving simulations, errors can be expected as will be shown later in this paper. All settings for traction controls and common model parameters are the same between the 2D and 3D locomotive models. One example of the simulated traction characteristics using the 2D model is also shown in Fig. 3. It can be seen that there are minor differences for the simulated adhesion characteristics between the 3D and 2D locomotive models. These minor differences are mainly generated from the DOF differences. In other words, the 3D model has more DOFs and more flexibility.

4 Simulation strategies

As the 2D and 3D locomotive models were developed in different software packages (in-house and GENSYS, respectively), different simulation strategies are used to combine them with the in-house LTD model for energy calculation studies. In this paper, the Message Passing Interface (MPI) technique [24] was used for the 2D locomotive model, whilst the Transmission Control Protocol/Internet Protocol (TCP/IP) technique [25] and the OpenMP technique [26] were used for the 3D locomotive model. Both MPI and OpenMP are able to conduct parallel computing to increase computing speeds.

4.1 For 2D locomotives

Both the LTD model and 2D locomotive model were developed in-house in the FORTRAN language. The communications between the two models and parallel computing were enabled by using the MPI technique. The simulation strategy is shown in Fig. 4. In simulations, multiple locomotives at arbitrary positions along the train can be considered. In each simulation step, the LTD model first computes traction references, locomotive speeds in LTD, coupler forces, resistance forces, and grade forces. In the LTD model, the traction references are traction forces before adhesion controls; the traction forces will be converted to traction torques and then go through adhesion controls. Traction torques after the adhesion controls will then be applied to individual wheelsets. Locomotive speeds in LTD are used as a synchronisation check between the two models. In other words, the simulated speeds of the locomotives in the two models should be the same. All this information from the LTD model is sent to corresponding locomotive models via the MPI communication protocol. Having received the information, MPI starts simulations of different locomotives simultaneously on multiple computer cores (parallel computing). In each locomotive model, external forces received from the LTD model are attached to corresponding bodies of the locomotive model. Vehicle system dynamics simulations with considerations of wheel–rail contact and traction control are conducted to determine the adhesion forces in wheel–rail interfaces. The sum of longitudinal adhesion forces from all wheel–rail contact is then calculated and sent to LTD as the traction force after adhesion control. This traction force after traction control is then attached to the corresponding locomotive models in the LTD model. LTD simulations then proceed to repeat the simulation cycle.

Fig. 4
figure 4

Simulation strategy for 2D locomotive and LTD models

Full size image

4.2 For 3D locomotives

The 3D locomotive model was developed using a commercial software package called GENSYS. A co-simulation interface is needed to facilitate the communications between the locomotive model and the LTD model. In this paper, a co-simulation interface based on TCP/IP was used. In addition, OpenMP was used to initiate multiple threads for parallel computing of the locomotive models, as shown in Fig. 5. The information exchanged between the LTD model and locomotive models is the same as those in Sect. 4.1 for the 2D locomotive models. Note that, in the 3D locomotive simulation strategy, TCP/IP is used as the communication technique as well as the synchronisation technique; the single function of OpenMP is to initiate multiple threads. Specifically, an interim module was developed between the LTD model and locomotive model. In each time step, the communication between the LTD model and the interim model is only handled by one single thread (Core 0); further actions will not be taken unless all information is sent, or all information is received. In this way, TCP/IP acts as a synchronisation function as well. The communications between the interim module and locomotive model are handled by multiple threads (Core 1–Core 3) that are initiated by OpenMP. In this way, the simulations of locomotive models can be conducted in the parallel computing fashion.

Fig. 5
figure 5

Simulation strategy for 3D locomotive and LTD models

Full size image

5 Energy calculation method and case studies

Having now achieved train dynamics simulations with considerations of wheel–rail contact and adhesion control, advanced energy calculations can be carried out. The energy consumption of each locomotive is simply calculated as

$$ E = \sum F_{\text{t}} S, $$
(4)

where \( E \) is the energy consumed; \( F_{\text{t}} \) is the traction force; and \( S \) is the translation distance during the time step. To compare energy calculation results using different models, three simulations were conducted:

  • Simulation 1: train dynamics simulation and energy calculation using 1D locomotive and LTD models without considerations of wheel–rail contact and adhesion control.

  • Simulation 2: train dynamics simulation and energy calculation using 2D locomotive and LTD models with considerations of the simplified wheel–rail contact model and adhesion control.

  • Simulation 3: train dynamics simulation and energy calculation using 3D locomotive and LTD models with considerations of the full wheel–rail contact model and adhesion control.

In these simulations, the simulated trains have the same train configuration with 3 locomotives + 105 wagons. The wagon models have mass properties of 120 t each and only one DOF each. The propulsion resistance formulation used in this paper is expressed as [19]

$$ F_{\text{pr}} = m_{\text{w}} \left( {2.943 + \frac{89.2}{{m_{\text{a}} }} + 0.0306V_{\text{km}} + \frac{{0.122V_{\text{km}}^{2} }}{{m_{\text{w}} }}} \right), $$
(5)

where \( m_{\text{w}} \) and \( m_{\text{a}} \) are vehicle total mass and axle load, respectively; and \( V_{\text{km}} \) is vehicle speed. Curving resistance formula is expressed as

$$ F_{\text{cr}} = m_{\text{w}} \left( {6116/R} \right), $$
(6)

where \( R \) is curve radius.

The simulated traction section is a flat section with a 300 m radius curve. The curve starts at 3.3 km and has the configuration of 50 m entry transition + 200 m circular curve + 50 m exit transition. Only traction was considered during the simulations; the trains start at 3.0 km with the speed of 5 km/h; locomotive traction starts at notch 1 and gradually increases to notch 8 with a 10-s interval. After this, traction notch was kept at notch 8 until the end of the simulations. Other notes for the simulations are as follows:

  • In all three simulations, a 0.35 adhesion limit was applied.

  • In Simulation 2 and Simulation 3, wheel slip limits were set as 7%.

  • In Simulation 2 and Simulation 3, adhesion model parameters were as listed in Table 1.

    Table 1 Adhesion model parameters
    Full size table

It is noted that the simulation cases in this paper focused on traction-only scenarios. These cases are used merely for demonstration purposes. Adhesion characteristics and studies for braking scenarios [27, 28] are also important. Wheel–rail contact and adhesion models used in this paper are also able to conduct braking studies.

Figure 6 presents the simulated traction forces from the three different models. Note that the horizontal axis of the figure is the location of the train on the track, i.e. the location of the leading locomotive of the train. As stated, the track section from 3.3 to 3.6 km is a section of curved track. The results show that there are traction force differences among different models in both tangent and curved sections. For example, at track location 3.15 km (tangent section), simulated traction forces from 1D, 2D, and 3D locomotives are 460, 400, and 416 kN, respectively. Using the result of the 3D locomotive model as a baseline, results from 1D and 2D locomotive models have 9.6% and − 3.8% differences, respectively. The result difference between 1D locomotive model and 3D locomotive model is that the static coefficient of friction is always larger than the coefficient of adhesion, as shown in Fig. 7. The coefficient of friction utilised by the locomotive, i.e. the coefficient of adhesion, is always lower than the static coefficient of friction. Regarding the relatively smaller difference between the 2D locomotive model and 3D locomotive model, this can be attributed to the lateral movements of the wheelsets as well as the overall model differences. The 3D locomotive model considers all DOFs of all major components of the simulated locomotive, whilst the 2D locomotive model does not consider lateral, roll, and yaw movements of the model components.

Fig. 6
figure 6

Simulated traction forces of the leading locomotive

Full size image
Fig. 7
figure 7

Illustrations of coefficient of friction and coefficient of adhesion

Full size image

The simulation results show that, at track location 3.45 km (curved section), simulated traction forces from 1D, 2D, and 3D locomotives are 460, 386, and 423 kN, respectively. Again, using the result of the 3D locomotive model as a baseline, results from 1D and 2D locomotive models have 8.7% and − 8.7% differences, respectively. It is noted that the difference between the 2D model and the 3D model is significantly larger than that in the tangent section. Deeper investigation into wheel–rail forces found that this larger difference resulted from the different number of contact points in these two simulations. Specifically, the 3D model had more than one wheel–rail contact point during curve negotiation, which generated extra traction forces. However, the 2D model only considers one contact point; therefore, the traction force did not increase during curve negotiation.

Figure 8 shows the simulated energy consumptions of the leading locomotive in three different simulations. The figure has also marked out the consumption values at 3.6 km (end of curved section) and 4.5 km (end of simulation) locations. It can be seen that evident differences are recorded during the simulations. Compared to the 3D model, differences in simulated energy consumptions of the 1D and 2D models within the 3.0–3.6 km section are 7.6% and − 4.9%, respectively; and the differences for the whole section (3.0–4.5 km) are 2.5% and − 1.7%, respectively.

Fig. 8
figure 8

Simulated energy consumptions of the leading locomotive

Full size image

Figure 9 shows the differences of energy consumption rates of different models at different locations. The consumption rates are directly related to locomotive traction forces at specific locations. In other words, a locomotive that has higher traction force also has higher energy consumption rates. The differences shown in Fig. 9 were calculated as the simulated traction force difference between a simpler model and the 3D model divided by the simulated traction force from the 3D model. Figure 9 shows that the major differences between the results of the 1D model and the 3D model were recorded at 3.2, 3.3, and 3.4 km where the locomotives have the maximum traction forces. The maximum difference was 12.5% at 3.3 km. The major differences between the results of the 2D model and the 3D model were recorded at 3.4 and 3.5 km where the locomotives were negotiating the curve. The maximum difference was − 8.6% at 3.4 km.

Fig. 9
figure 9

Differences of simulated energy consumption rates

Full size image

6 Conclusions

Traction force simulation for locomotives is the prerequisite of train energy studies; it can be regarded as a fundamental part of such research topic. Traction forces are known to be influenced by wheel–rail adhesions, locomotive traction control, and locomotive dynamics; however, none of these factors can be considered by conventional train energy calculation models that simplify locomotives as single DOF rigid bodies or as part of a mass block that is used to simulate the train as a whole.

Two train models that have different levels of complexities were developed in this paper to fill the stated research gap. The first model, which is the simpler one, uses a 2D locomotive model that has 27 DOFs and a simplified wheel–rail contact model. The second model uses a 3D locomotive model with 54 DOFs and a fully detailed wheel–rail contact model. Both models were integrated into an LTD model with the consideration of locomotive adhesion control. Different simulation strategies were used to cater for the differences in programming languages and software packages. Specifically, the 2D locomotives were simulated with the LTD model by using the MPI technique; the 3D counterparts were simulated by using OpenMP and TCP/IP techniques.

Comparative studies were conducted to simulate energy consumptions of a heavy haul train by using a conventional model used in LTD simulations and these two new models. These three models are called 1D model, 2D model, and 3D model, respectively, in the paper. Simulation results show that the 3D model reports less energy consumption than the 1D model due to the consideration of the wheel–rail adhesion model and traction control in the 3D model. The maximum difference in energy consumption rate between the 3D model and the 1D model was 12.5%. Due to the consideration of multiple contact points in the 3D model, it reports higher energy consumption than the 2D model. An 8.6% maximum difference in energy consumption rate was reported during curve negotiation between the 3D model and the 1D model. The results also show that the differences in energy calculations are direct reflections of the differences of the simulated traction forces. This confirms that, to achieve more accurate energy consumption calculations, more accurate traction force simulations are required.