Case study specification#
Written by Marc Budinger, INSA Toulouse, France
We will now develop tools allowing:
generate mission profiles corresponding to a new need for an urban transport line
optimize certain mission profiles choices to minimize consumption
In order to have flexibility in the definition of the case studies, an object-oriented programming approach will be used.
Electric vehicles#
Different types of urban electric vehicles exist (tram, trolley bus, electric bus), we can define a general vehicle class wich can represent the multiple possible supports.
class vehicle:
def __init__(self, MaxSpeed, MaxAcc,NumberPass,Weight,Crr,Cd,FrontArea):
self.MaxSpeed = MaxSpeed # [m/s] maximal speed
self.MaxAcc = MaxAcc # [m/s²] maximal acceleration
self.NumberPass = NumberPass # [-] passenger capacity (3p/m²)
self.Weight = Weight # [kg] effective weight (with passenger 62kg/p)
self.Crr = Crr # [-] rolling resistance coefficient
self.Cd = Cd # [-] drag coefficient
self.FrontArea = FrontArea # [m²] frontal area
self.Fmax = Weight*MaxAcc # [N] max traction force
self.Pmax = self.Fmax*self.MaxSpeed # [W] max corner power
Here you will find some useful vehicle characteristics for the rest of the study.
Attribute |
Tram |
Trolleybus |
Bus |
Car (Tesla 3) |
|---|---|---|---|---|
Maximum operational speed |
60 km/h |
60 km/h |
50 km/h |
100 km/h (1st gear) |
Acceleration and braking |
1.2 m/s² |
1.2 m/s² |
1.2 m/s² |
4.6 m/s² |
Passenger capacity @ 3p/m² |
220 |
138 |
95 |
4 |
Vehicle effective weight @ 62 kg/p |
57049 kg |
27656 kg |
19500 kg |
1800 kg |
Rolling resistance coefficient |
0.006 |
0.015 |
0.015 |
0.015 |
Aerodynamic drag coefficient |
0.6 |
0.6 |
0.6 |
0.23 |
Frontal area |
8.5 m² |
8.5 m² |
8.5 m² |
2.22 m² |
Exercice: Define 4 instances of the class
vehiclecorresponding to the vehicles above.
Tram = vehicle(MaxSpeed = 60*1e3/3600, MaxAcc=1.2,NumberPass=220,Weight=57049,Crr=0.006,Cd=0.6,FrontArea=8.5)
TrolleyBus = vehicle(60*1e3/3600,1.2,138,27656,0.015,0.6,8.5)
Bus = vehicle(50*1e3/3600,1.2,95,19500,0.015,0.6,8.5)
Car = vehicle(60*1e3/3600,4.6,4,1800,0.015,0.23,2.22)
Generation and mission profile simulation#
The objective now is to simulate the dynamic evolution of the main mechanical quantities (position, speed, acceleration, tractive effort, driving power) on typical sections of urban routes.
When a vehicle travels from point A (Starting point) to point B (Stopping point ) it generally passes through four stages namely:
Accelerating mode
Constant speed mode
Coasting/freewheel mode
Braking mode
The control points for switching from one mode to another are defined in the following code using 3 parameters \(k_1, k_2, k_3\) as shown in the figure below.
The following code implements a control logic comprising the 4 phases defined previously and a numerical integration of the differential equation corresponding to the longitudinal dynamics of the vehicle.
Question: What evolution should be added to take into account variations in altitude?
import numpy as np
from scipy import signal
from scipy.optimize import fmin_slsqp
from scipy.optimize import differential_evolution
import matplotlib.pyplot as plt
from scipy.integrate import odeint
class SimulSection:
def __init__(self,Vehicle,Distance, MeanSpeed,BrakeRatioMax,dt):
# Parameters definition
self.Vehicle = Vehicle # vehicle parameters
self.Distance = Distance # [m] distance to travel
self.MeanSpeed = MeanSpeed # [m/s] Mean Speed < MaxSpeep
self.TravelTime = self.Distance / self.MeanSpeed # [s] Travel time
self.BrakeRatioMax = BrakeRatioMax # [-], inferieur à 1, ratio de puissance de freinage / Puisssance d'acceleration
self.rho = 1.25 # [kg/m3] air density
self.g = 9.81 # [m/s²] gravity
self.dt = dt # Time step for numerical solver
# Tests sur les lois de mouvement
# https://fr.wikipedia.org/wiki/Loi_de_mouvement
if (MeanSpeed>=Vehicle.MaxSpeed):
print("Mean Speed is to high : Mean Speed > Vehicle Max Speed")
# profil triangulaire : minimise l'acceleration
# amax = 4⋅xf/T2 ; vmax = 2⋅xf/T ; amax = vmax / (T/2) ; vmean = vmax / 2
# amax = 4*vmean^2 / xf
# Calcul de l'acceleration possible a pleine vitesse
AmaxFullSpeed=(self.Vehicle.Fmax - 1/2*self.rho*self.Vehicle.Cd*self.Vehicle.FrontArea*self.Vehicle.MaxSpeed**2 - self.Vehicle.Weight*self.g*self.Vehicle.Crr)/self.Vehicle.Weight
if (AmaxFullSpeed < 4*self.MeanSpeed**2 / self.TravelTime):
print("Vehicle Max Acceleration (depending of max force) is too low or travel distance too small or mean speed too small")
# dynamic model for acceleration, coasting and braking
def model(self, y,F):
# state
x, dxdt = y
# System: acceleration calculation
dxdt2 = (F-self.Vehicle.Crr*self.Vehicle.Weight*self.g*np.sign(dxdt)
-1/2*self.rho*self.Vehicle.Cd*self.Vehicle.FrontArea*dxdt**2**np.sign(dxdt))/self.Vehicle.Weight
dydt = [dxdt, dxdt2]
return dydt
# https://perso.crans.org/besson/publis/notebooks/Runge-Kutta_methods_for_ODE_integration_in_Python.html
# solver numerique
def solver(self, x):
k1, k2, k3 = x
dt=self.dt # [s] pas de temps pour l'integration
t=0
self.tsection = [0]
self.xsection = [0]
self.vsection = [0]
self.asection = [0]
self.psection = [0]
self.Fsection = [0]
self.NRJsection = [0]
NRJ=0
y0 = np.array([0, 0])
y= y0
while (y[1]>=0):
# Traction/Braking force
if (y[0]<k1*k2*k3*self.Distance):
F=self.Vehicle.Fmax
dydt = self.model(y, F)
elif (y[0]<k3*k2*self.Distance):
F=self.Vehicle.Crr*self.Vehicle.Weight*self.g+1/2*self.rho*self.Vehicle.Cd*self.Vehicle.FrontArea*y[1]**2
dydt = [y[1], 0]
elif (y[0]<k3*self.Distance):
F=0
dydt = self.model(y, F)
else:
F=-self.BrakeRatioMax*self.Vehicle.Fmax
dydt = self.model(y, F)
self.Fsection = self.Fsection + [F]
# Euler integration de la postion y[0] et de la vitesse y[1]
y = y + dt * np.array(dydt)
t = t + dt
NRJ = NRJ+y[1]*F*dt
self.tsection = self.tsection + [t]
self.xsection = self.xsection + [y[0]]
self.vsection = self.vsection + [y[1]]
self.psection = self.psection + [y[1]*F]
self.asection = self.asection + [dydt[1]]
self.NRJsection = self.NRJsection + [NRJ]
def plot(self):
fig, axs = plt.subplots(3,2)
axs[0,0].plot(self.tsection,self.xsection,'b-',label='Simulation')
axs[0,0].plot(self.tsection,self.Distance*np.ones(len(self.tsection)),'g-',label='Specification')
axs[0,0].set_ylabel("Position (m)")
axs[0,0].legend()
axs[0,0].grid()
axs[1,0].plot(self.tsection,self.vsection,'b-', label='Simulation')
axs[1,0].plot(self.tsection,self.MeanSpeed*np.ones(len(self.tsection)),'g-',label='Specification')
axs[1,0].plot(self.tsection,np.mean(self.vsection)*np.ones(len(self.tsection)),'r--',label='Mean')
axs[1,0].set_ylabel("Speed (m/s)")
#axs[1].sharex(axs[0])
axs[1,0].grid()
axs[1,0].set_xlabel('Time (s)')
axs[0,1].plot(self.tsection,self.Fsection,'b-')
axs[0,1].set_ylabel("Force (N)")
axs[0,1].grid()
axs[2,0].plot(self.tsection,self.asection,'b-')
axs[2,0].set_ylabel("Acceleration (m/s²)")
axs[2,0].grid()
axs[2,0].set_xlabel('Time (s)')
axs[1,1].plot(self.tsection,np.array(self.psection)*1e-3,'b-')
axs[1,1].set_ylabel("Power (kW)")
axs[1,1].grid()
axs[1,1].set_xlabel('Time (s)')
axs[2,1].plot(self.tsection,np.array(self.NRJsection)*1e-3,'b-')
axs[2,1].set_ylabel("Energy (kJ)")
axs[2,1].grid()
axs[2,1].set_xlabel('Time (s)')
fig.tight_layout()
def ConsumptionPerPax(self,x):
self.solver(x)
print("Consumption per passenger: %.2f kJ/(Pax.km)"%(self.NRJsection[-1]/self.Vehicle.NumberPass/self.Distance))
return (self.NRJsection[-1]/self.Vehicle.NumberPass/self.Distance)
def MaxEnergyPoint(self,x):
self.solver(x)
print("Max energy discharge: %.0f kJ"%(max(self.NRJsection)*1e-3))
return max(self.NRJsection)*1e-3
Example of use of the previous code:
# Definition of the section under study
Trajet_Tram = SimulSection(Vehicle=Tram, Distance=670, MeanSpeed=26*1e3/3600,BrakeRatioMax=1.0,dt=1)
Trajet_Car = SimulSection(Vehicle=Car, Distance=670, MeanSpeed=26*1e3/3600,BrakeRatioMax=1.0,dt=1)
print('Duration of traject for given mean speed: ' + str(670/(26*1e3/3600)) +'s')
# Profile integration, plotting and consummption per passenger evaluation
X=[0.05, 1, 0.95]
Trajet_Tram.solver(X)
Trajet_Tram.plot()
Trajet_Tram.ConsumptionPerPax(X)
Trajet_Tram.MaxEnergyPoint(X)
# Profile integration, plotting and consummption per passenger evaluation
X=[0.005, 1, 0.999]
Trajet_Car.solver(X)
Trajet_Car.plot()
Trajet_Car.ConsumptionPerPax(X)
Vehicle Max Acceleration (depending of max force) is too low or travel distance too small or mean speed too small
Duration of traject for given mean speed: 92.76923076923077s
Consumption per passenger: 21.18 kJ/(Pax.km)
Max energy discharge: 5039 kJ
Consumption per passenger: 103.06 kJ/(Pax.km)
103.06442056207767
Two extrema cases are :
a triangular profile (acceleration/deceleration, without constant speed & freewheel)
a trapezoidal profile (without freewheel)
Question 1: Give the relationship between acceleration, average speed and distance traveled for triangular profile if we assume the acceleration and braking phases are perfectly linear and of the same slope (i.e. resistive force negligible): does it match a warning message of the SimulSection class?
Question 2: Give the relationship between k1 (k2=1), average speed and distance traveled if we assume the acceleration and braking phases are perfectly linear and of the same slope.
# Calculus of lower acceleration
a_min = 4*(26*1e3/3600)**2/670
print('Min acceleration is: ' + str(a_min) +'m/s²')
Tram_low_acc = vehicle(MaxSpeed = 60*1e3/3600, MaxAcc=a_min,NumberPass=220,Weight=57049,Crr=0.006,Cd=0.6,FrontArea=8.5)
# Triangular profile test
X=[1, 1, 0.6] #if resistive force not negligieble in front of propulsive force not k1>0.5
Trajet_Tram_low_acc = SimulSection(Vehicle=Tram_low_acc, Distance=670, MeanSpeed=26*1e3/3600,BrakeRatioMax=1.0,dt=1)
Trajet_Tram_low_acc.solver(X)
Trajet_Tram_low_acc.plot()
Trajet_Tram_low_acc.ConsumptionPerPax(X)
Trajet_Tram_low_acc.MaxEnergyPoint(X)
Min acceleration is: 0.3114059332964806m/s²
Vehicle Max Acceleration (depending of max force) is too low or travel distance too small or mean speed too small
Consumption per passenger: 18.74 kJ/(Pax.km)
Max energy discharge: 7420 kJ
7419.990488324914
# Parameters calculation
# on pose x = sqrt(k1)
from scipy.optimize import root
a = Trajet_Tram.MeanSpeed/(2*Trajet_Tram.Vehicle.MaxAcc*Trajet_Tram.Distance)**0.5
f = lambda x: a+(x*(1-x))**1.5-(x*(1-x))**0.5
X=np.linspace(0,0.5,100)
plt.plot(X,f(X),'b-')
plt.grid()
plt.show()
sol=root(f, 0.1)
k1=float(sol.x)
print(k1)
0.036125972011676784
/tmp/ipykernel_2603/2473967107.py:17: DeprecationWarning: Conversion of an array with ndim > 0 to a scalar is deprecated, and will error in future. Ensure you extract a single element from your array before performing this operation. (Deprecated NumPy 1.25.)
k1=float(sol.x)
k3=1-k1
X=[k1,1,k3]
print(X)
Trajet_Tram.solver(X)
Trajet_Tram.plot()
Trajet_Tram.ConsumptionPerPax(X)
Trajet_Tram.MaxEnergyPoint(X)
[0.036125972011676784, 1, 0.9638740279883232]
Consumption per passenger: 20.44 kJ/(Pax.km)
Max energy discharge: 4450 kJ
4449.63595287904
Questions: Is this type of profile the one that minimizes the energy consumed ? What rules do you recommend to minimize the energy consumed?
Optimizing power consumption using speed profile#
In a more general case and contrary to the previous hypothesis, the braking power is not necessarily identical to that possible during acceleration. The following article shows that it is possible to optimize the consumption of a tram:
Tian, Z., Zhao, N., Hillmansen, S., Roberts, C., Dowens, T., & Kerr, C. (2019). SmartDrive: Traction energy optimization and applications in rail systems. IEEE Transactions on Intelligent Transportation Systems, 20(7), 2764-2773.
The previous class is now extended to allow optimization of variables \(k_1, k_2, k_3\) in order to minimize energy consumption (even with braking power lower than acceleration power).
class OptimSection(SimulSection):
# Gradient optimization
def objectifConso(self, x):
self.solver(x)
NRJmin = self.xsection[-1]*self.Vehicle.Weight*self.g*self.Vehicle.Crr
return (self.NRJsection[-1]/NRJmin)
def contraintes(self, x):
self.solver(x)
TimeTrajet = self.Distance / self.MeanSpeed
PowerAccBrakeRatio = abs(min(self.psection))/max(self.psection)
print('Brake power/Acceleration power ratio is: '+ str(PowerAccBrakeRatio))
# contraintes en distance, vitesse moyenne (ou duree de trajet, vitesse au freinage), ratio puissance freinage / acc
return [(self.xsection[-1]-self.Distance)/self.Distance, (TimeTrajet-self.tsection[-1])/TimeTrajet,
self.BrakeRatioMax-PowerAccBrakeRatio ]
def optimizeConso(self,x0):
Xbound= [(0.0, 0.3), (0.1, 1), (0.5,.99)]
Xopt=fmin_slsqp(func=self.objectifConso, x0=x0, f_ieqcons = self.contraintes, bounds=Xbound, epsilon=1e-2)
return Xopt
# Differential evolution optimization
def objectifGConso(self, x):
self.solver(x)
NRJmin = self.xsection[-1]*self.Vehicle.Weight*self.g*self.Vehicle.Crr
pen=0
TimeTrajet = self.Distance / self.MeanSpeed
PowerAccBrakeRatio = abs(min(self.psection))/max(self.psection)
VecC=[(self.xsection[-1]-self.Distance)/self.Distance, (TimeTrajet-self.tsection[-1])/TimeTrajet, self.BrakeRatioMax-PowerAccBrakeRatio ]
for C in VecC:
if (C<0):
pen=pen-1e2*C
return (self.NRJsection[-1]/NRJmin+pen)
def optimizeGConso(self,x0):
Xbound= [(0.0, 0.3), (0.1, 1), (0.5,0.99)]
res=differential_evolution(func=self.objectifGConso, x0=x0, bounds=Xbound)
print(res)
return res.x
# Differential evolution optimization
def objectifGEmax(self, x):
self.solver(x)
NRJmin = self.xsection[-1]*self.Vehicle.Weight*self.g*self.Vehicle.Crr
pen=0
TimeTrajet = self.Distance / self.MeanSpeed
PowerAccBrakeRatio = abs(min(self.psection))/max(self.psection)
VecC=[(self.xsection[-1]-self.Distance)/self.Distance, (TimeTrajet-self.tsection[-1])/TimeTrajet, self.BrakeRatioMax-PowerAccBrakeRatio ]
for C in VecC:
if (C<0):
pen=pen-1e2*C
return (max(self.NRJsection)/NRJmin+pen)
def optimizeGEmax(self,x0):
Xbound= [(0.0, 0.3), (0.1, 1), (0.5,0.99)]
res=differential_evolution(func=self.objectifGEmax, x0=x0, bounds=Xbound)
print(res)
return res.x
The usage of the class is quite simular that the previous one, except an optimizewhich enable to optimize the \(k_1, k_2, k_3\) set of profile parameter:
# Definition of the section under study
Trajet = OptimSection(Vehicle=Tram, Distance=670, MeanSpeed=26*1e3/3600,BrakeRatioMax=0.6,dt=0.25)
# Initial variables vector
X=[0.036125972011676784, 1, 0.9638740279883232]
# Optimization of the consumption
Xopt1=Trajet.optimizeGConso(X) # Global optimization with differential evolution algorithm
#Xopt1=Trajet.optimizeConso(X) # Opitmization with gradient algorithm
# print and plot results
print("Optimal vector :", Xopt1)
print("Constraints vector :", Trajet.contraintes(Xopt1))
Trajet.solver(Xopt1)
Trajet.plot()
Trajet.ConsumptionPerPax(Xopt1)
Trajet.MaxEnergyPoint(Xopt1)
# Optimization of the max energy point
print('\n######################################\n')
Xopt2=Trajet.optimizeGEmax(X) # Global optimization with differential evolution algorithm
# print and plot results
print("Optimal vector :", Xopt2)
print("Constraints vector :", Trajet.contraintes(Xopt2))
Trajet.solver(Xopt2)
Trajet.plot()
Trajet.ConsumptionPerPax(Xopt2)
Trajet.MaxEnergyPoint(Xopt2)
Vehicle Max Acceleration (depending of max force) is too low or travel distance too small or mean speed too small
message: Optimization terminated successfully.
success: True
fun: 1.105402086832088
x: [ 6.035e-02 7.181e-01 9.644e-01]
nit: 22
nfev: 1039
population: [[ 6.035e-02 7.181e-01 9.644e-01]
[ 4.553e-02 9.683e-01 9.422e-01]
...
[ 5.829e-02 7.570e-01 9.638e-01]
[ 5.535e-02 7.730e-01 9.547e-01]]
population_energies: [ 1.105e+00 1.111e+00 ... 1.106e+00 1.108e+00]
Optimal vector : [0.06034734 0.71808202 0.96436532]
Brake power/Acceleration power ratio is: 0.4779240994786292
Constraints vector : [0.01063718582405839, 0.00020729684908794098, 0.1220759005213708]
Consumption per passenger: 17.05 kJ/(Pax.km)
Max energy discharge: 3692 kJ
######################################
message: Optimization terminated successfully.
success: True
fun: 1.4627161460359577
x: [ 2.666e-01 2.197e-01 9.695e-01]
nit: 46
nfev: 2119
population: [[ 2.666e-01 2.197e-01 9.695e-01]
[ 2.517e-01 2.352e-01 9.696e-01]
...
[ 2.692e-01 2.199e-01 9.697e-01]
[ 2.672e-01 2.206e-01 9.685e-01]]
population_energies: [ 1.463e+00 1.474e+00 ... 1.463e+00 1.463e+00]
Optimal vector : [0.26655047 0.21971369 0.96946951]
Brake power/Acceleration power ratio is: 0.33307979392645
Constraints vector : [0.0006461216493124681, 0.00020729684908794098, 0.26692020607355]
Consumption per passenger: 16.97 kJ/(Pax.km)
Max energy discharge: 3293 kJ
3292.9388688877516
Question: Compare this optimized profile to a profile with only 3 segments: maximum acceleration, constant speed, maximum deceleration. Evaluate effect of specifications (distance, braking power/acceleration power ratio).
Different vehicules comparison#
Exercice: For the same travel (distance, main speed) compare the consumption of different vehicles (tramway, tram, bus, car). The function
ConsumptionPerPaxreturn and print the consumption per km and per passenger. Compare the results with this publication.
print("Calculation of energy consumption of different vehicles: ")
print("")
print("Tramway :")
Trajet = OptimSection(Tram, 670, 26*1e3/3600,0.5,1)
X=[0.036125972011676784, 1, 0.9638740279883232]
Xopt=Trajet.optimizeGConso(X)
Trajet.ConsumptionPerPax(Xopt)
print("----")
print("Trolleybus :")
Trajet = OptimSection(TrolleyBus, 670, 26*1e3/3600,0.5,1)
X=[0.036125972011676784, 1, 0.9638740279883232]
Xopt=Trajet.optimizeGConso(X)
Trajet.ConsumptionPerPax(Xopt)
print("----")
print("Bus :")
Trajet = OptimSection(Bus, 670, 26*1e3/3600,0.5,1)
X=[0.036125972011676784, 1, 0.9638740279883232]
Xopt=Trajet.optimizeGConso(X)
Trajet.ConsumptionPerPax(Xopt)
print("----")
print("Car :")
Trajet = OptimSection(Car, 670, 26*1e3/3600,0.5,1)
X=[0.036125972011676784, 1, 0.9638740279883232]
Xopt=Trajet.optimizeGConso(X)
Trajet.plot()
Trajet.ConsumptionPerPax(Xopt)
Calculation of energy consumption of different vehicles:
Tramway :
Vehicle Max Acceleration (depending of max force) is too low or travel distance too small or mean speed too small
message: Optimization terminated successfully.
success: True
fun: 1.2482194479805129
x: [ 8.460e-02 4.903e-01 9.752e-01]
nit: 20
nfev: 949
population: [[ 8.460e-02 4.903e-01 9.752e-01]
[ 8.721e-02 4.809e-01 9.679e-01]
...
[ 8.340e-02 4.973e-01 9.721e-01]
[ 8.307e-02 4.870e-01 9.724e-01]]
population_energies: [ 1.248e+00 1.249e+00 ... 1.248e+00 1.248e+00]
Consumption per passenger: 19.71 kJ/(Pax.km)
----
Trolleybus :
Vehicle Max Acceleration (depending of max force) is too low or travel distance too small or mean speed too small
message: Optimization terminated successfully.
success: True
fun: 1.119398008169179
x: [ 5.291e-02 8.126e-01 9.812e-01]
nit: 20
nfev: 949
population: [[ 5.291e-02 8.126e-01 9.812e-01]
[ 5.494e-02 8.007e-01 9.779e-01]
...
[ 5.547e-02 8.013e-01 9.767e-01]
[ 1.214e-01 6.640e-01 9.329e-01]]
population_energies: [ 1.119e+00 1.120e+00 ... 1.120e+00 1.186e+00]
Consumption per passenger: 33.80 kJ/(Pax.km)
----
Bus :
Vehicle Max Acceleration (depending of max force) is too low or travel distance too small or mean speed too small
message: Optimization terminated successfully.
success: True
fun: 1.1382110888204808
x: [ 5.807e-02 7.584e-01 9.701e-01]
nit: 19
nfev: 904
population: [[ 5.807e-02 7.584e-01 9.701e-01]
[ 4.432e-02 9.014e-01 9.362e-01]
...
[ 4.677e-02 8.439e-01 9.638e-01]
[ 4.648e-02 8.073e-01 9.697e-01]]
population_energies: [ 1.138e+00 1.149e+00 ... 1.145e+00 1.142e+00]
Consumption per passenger: 34.50 kJ/(Pax.km)
----
Car :
message: Optimization terminated successfully.
success: True
fun: 1.3512967382800352
x: [ 8.441e-03 7.042e-01 9.861e-01]
nit: 38
nfev: 1759
population: [[ 8.441e-03 7.042e-01 9.861e-01]
[ 8.516e-03 7.052e-01 9.860e-01]
...
[ 5.646e-03 7.039e-01 9.857e-01]
[ 8.027e-03 6.981e-01 9.877e-01]]
population_energies: [ 1.351e+00 1.351e+00 ... 1.351e+00 1.351e+00]
Consumption per passenger: 89.63 kJ/(Pax.km)
89.63150769927998
Homework#
Read the following paper:
Jaafar, A., Sareni, B., Roboam, X., & Thiounn-Guermeur, M. (2010, September). Sizing of a hybrid locomotive based on accumulators and ultracapacitors. In 2010 IEEE Vehicle Power and Propulsion Conference (pp. 1-6). IEEE. [pdf]