Virtual Simulation of Electric Bus Fleets for City Bus Transport Electriﬁcation Planning

: City bus transport electriﬁcation has a strong potential of improving city air quality, reducing noise pollution and increasing passenger satisfaction. Since the city bus operation is rather deterministic and intermittent, the driving range- and charging-related concerns may be e ﬀ ectively overcome by means of fast charging at end stations and / or slow charging in depot. In order to support decision making processes, a simulation tool for planning of city bus transport electriﬁcation has been developed and it is presented in this paper. The tool is designed to use real / recorded driving cycles and techno-economic data, in order to calculate the optimal type and number of e-buses and chargers, and predict the total cost of ownership including investment and exploitation cost. The paper focuses on computationally e ﬃ cient e-bus ﬂeet simulation including powertrain control and charging management aspects, which is illustrated through main results of a pilot study of bus transport electriﬁcation planning for the city of Dubrovnik.


Introduction
Due to environmental concerns, there is a strong tendency of electrifying road transport systems by means of introducing different types of electric vehicles [1]. Apart from reducing pollutant and CO 2 emissions, electric vehicles (EV) are characterised by substantially reduced noise pollution, lower operating cost (including energy and maintenance cost) and generally better driving characteristics. On the other hand, higher investment cost, slow battery charging and limited driving range inhibit their faster proliferation [2]. This is why the transition to fully electric vehicles (FEV or BEV) is characterised by application of hybrid electric vehicles (HEV) and plug-in hybrid electric vehicles (PHEV) [3].
City bus transport is a natural candidate for electrification, aiming at improvement in city air quality and reduction of noise. Since the city bus routes are known in advance and the operation is intermittent, the range-and charging-related issues are of lesser significance than with passenger cars. These issues can be tackled by the following two basic approaches [4]: (i) the buses are equipped with large enough battery packs to sustain half a day or even full day of operation, and the buses are efficiently recharged by using slow charging (typically overnight); and (ii) the battery size is minimised and superfast charging is employed at bus stops (typically end stations). Therefore, it is generally of interest to find optimal locations of charging stations, as considered in [5][6][7][8][9][10][11] with a focus on passenger cars and urban areas.
In [12], e-buses and corresponding charging systems are analysed from the total cost perspective by using data related to route specifications, timetables and other local conditions. Additionally, the authors have developed a user-friendly tool which enables the user to investigate and quantify Figure 2 illustrates the structure of the developed tool, including the connections between the main modules. The tool is driven by recorded driving cycle data, and as the main output, it delivers the TCO over the projected fleet operational period (12 years, herein). The DPPM transforms the recorded driving data into individual driving cycles, and it also calculates various statistical features characterising the conventional city bus transport behaviours. The EBSM provides computer simulations of different types of city buses (CONV, HEV, PHEV, BEV) over the driving cycles extracted by the DPPM. The module outputs include the individual bus energy consumption (fuel and/or electricity) and various features of powertrain response (e.g., engine/e-motor operating points, gear ratio trajectories, etc.). The COM utilises the outputs of DPPM and EBSM to simulate the overall city bus fleet over the recorded driving cycles and optimise the PHEV-and BEV-type bus charging configuration and management. This module provides the number, location and type of chargers, the bus battery capacity and the number of reserve buses in the BEV case, which are required to fulfil the driving routes with sufficient battery charge. The COM also outputs the total fuel and/or electricity consumption over the considered period of operation. The TEAM uses the output data from the COM module, as well as the data on bus transport investment and exploitation/maintenance cost, in order to calculate the TCO.
The simulation tool is written in Python object-oriented programming language, with computationally demanded routines coded in C language. The Python language has been selected because it is available in open-source domain, it supports multiple platforms, has a vast number of available libraries, and is dynamically typed language (i.e., does not require declaring of data type of a variable) convenient for rapid prototyping. The simulation tool is designed in a user-friendly way (based on a graphical user interface (GUI) including windows, tabs, I/O data interfaces, etc.) and has    Figure 2 illustrates the structure of the developed tool, including the connections between the main modules. The tool is driven by recorded driving cycle data, and as the main output, it delivers the TCO over the projected fleet operational period (12 years, herein). The DPPM transforms the recorded driving data into individual driving cycles, and it also calculates various statistical features characterising the conventional city bus transport behaviours. The EBSM provides computer simulations of different types of city buses (CONV, HEV, PHEV, BEV) over the driving cycles extracted by the DPPM. The module outputs include the individual bus energy consumption (fuel and/or electricity) and various features of powertrain response (e.g., engine/e-motor operating points, gear ratio trajectories, etc.). The COM utilises the outputs of DPPM and EBSM to simulate the overall city bus fleet over the recorded driving cycles and optimise the PHEV-and BEV-type bus charging configuration and management. This module provides the number, location and type of chargers, the bus battery capacity and the number of reserve buses in the BEV case, which are required to fulfil the driving routes with sufficient battery charge. The COM also outputs the total fuel and/or electricity consumption over the considered period of operation. The TEAM uses the output data from the COM module, as well as the data on bus transport investment and exploitation/maintenance cost, in order to calculate the TCO.

Organisational Structure of Simulation Tool
Energies 2020, 13, 3410 4 of 23 main storage for recorded driving cycle data, and plays the role of an intermediary between the main tool modules. In addition, the simulation tool includes the Data Management Module (DMM), which provides greater flexibility and adaptability to different cities' transport system configurations. DMM enables the user to define all static data (system parameters) required by the simulation tool, e.g., those related to vehicle model parameters, end-station and depot locations, charging station parameters and techno-economic data.

General Description
According to the DPPM flowchart shown in Figure 3, the recorded driving data are first processed and stored in a database. Next, the driving cycles are extracted by using the DMM-based data related to geographical coordinates of end stations and depot (Subsection 3.2). At the same time, the fleet statistics are calculated for the entire fleet and individual buses (Subsection 3.3). Finally, the module outputs including driving cycles and statistical features can be plotted in different formats, and they are saved into the database.  The simulation tool is written in Python object-oriented programming language, with computationally demanded routines coded in C language. The Python language has been selected because it is available in open-source domain, it supports multiple platforms, has a vast number of available libraries, and is dynamically typed language (i.e., does not require declaring of data type of a variable) convenient for rapid prototyping. The simulation tool is designed in a user-friendly way (based on a graphical user interface (GUI) including windows, tabs, I/O data interfaces, etc.) and has in mind the transferability to other cities using a common/shared database. The database serves as a main storage for recorded driving cycle data, and plays the role of an intermediary between the main tool modules. In addition, the simulation tool includes the Data Management Module (DMM), which provides greater flexibility and adaptability to different cities' transport system configurations. DMM enables the user to define all static data (system parameters) required by the simulation tool, e.g., those related to vehicle model parameters, end-station and depot locations, charging station parameters and techno-economic data.  The driving cycle segmentation process resulted in a total of 122,727 extracted driving cycles. An example of recorded driving cycle is shown in Figure 4. Note that the vehicle velocity is directly measured, while the road grade is reconstructed from the rate of change of altitude when expressed with respect to the travelled distance. In order to reduce the noise in the reconstructed road slope profile, the altitude signal is pre-filtered by a low-pass double-sided Butterworth filter [16].

Extraction of Driving Cycles
A single driving cycle is defined by the velocity vs. time and road slope vs. travelled distance profiles between two consecutive end stations, including depot (see Figure 1). The corresponding time profiles of cumulative fuel consumption are also extracted.
The driving cycle segmentation process resulted in a total of 122,727 extracted driving cycles. An example of recorded driving cycle is shown in Figure 4. Note that the vehicle velocity is directly measured, while the road grade is reconstructed from the rate of change of altitude when expressed with respect to the travelled distance. In order to reduce the noise in the reconstructed road slope profile, the altitude signal is pre-filtered by a low-pass double-sided Butterworth filter [16].

Calculation of Vehicle Fleet Statistics
A rich set of statistically significant driving features is calculated for the purpose of actual/conventional city-bus transport system characterisation and in support of transport electrification (e.g., locating charging stations, Section 5). The features related to individual buses, all given per-day-basis, include: the total fuel consumption and distance travelled; average fuel consumption in L/100 km; the total time the bus is dwelling at depot, individual end stations or any other locations (typically bus stops); total driving time; mean velocity; number of bus stops per kilometre; number of bus visits to depot and end stations. The results are stored in a two-dimensional (2D) matrix (one per bus), whose rows and columns represent individual days and the statistical features, respectively.
Once the individual statistics are stored, they are further used to calculate the same features for the entire fleet on the basis of individual day, week, month or year. The results can be presented in different ways, e.g., instead of individual dwelling times, one can get the information about percentage share of time the buses or entire fleet are resting at depot, end stations and other locations (see Figure 5). Specific fleet features requiring denser sampling (30 min, herein) are also calculated.

Calculation of Vehicle Fleet Statistics
A rich set of statistically significant driving features is calculated for the purpose of actual/conventional city-bus transport system characterisation and in support of transport electrification (e.g., locating charging stations, Section 5). The features related to individual buses, all given per-day-basis, include: the total fuel consumption and distance travelled; average fuel consumption in L/100 km; the total time the bus is dwelling at depot, individual end stations or any other locations (typically bus stops); total driving time; mean velocity; number of bus stops per kilometre; number of bus visits to depot and end stations. The results are stored in a two-dimensional (2D) matrix (one per bus), whose rows and columns represent individual days and the statistical features, respectively.
Once the individual statistics are stored, they are further used to calculate the same features for the entire fleet on the basis of individual day, week, month or year. The results can be presented in Energies 2020, 13, 3410 6 of 23 different ways, e.g., instead of individual dwelling times, one can get the information about percentage share of time the buses or entire fleet are resting at depot, end stations and other locations (see Figure 5). Specific fleet features requiring denser sampling (30 min, herein) are also calculated. An example of the daily average fleet velocity profile is shown in Figure 6. Other features related to entire fleet include: average number of buses being parked at depot on 24-hour time basis (Figure 6), clusters of buses parking durations in relation to geographical coordinates, count of transitions between individual end stations/depot, etc.
Energies 2020, 13, x FOR PEER REVIEW 6 of 25 rarely visit the depot in other time intervals, particularly over the work days (when the average number of buses is lower than 1).   Figure 5 confirms that the share of total time of buses being parked in the depot is relatively small (approx. 30%) and comparable to the share of end-station parking time (approx. 25%). The rest of the time the buses spend in driving (40%), while only for a small portion of time (approx. 5%), they rest elsewhere, typically at bus stops. Figure 7a indicates that there are significant differences in bus resting time at different end stations. The average resting durations for most pronounced end stations are between 10 and 20 min, thus making them good candidates for installation of fast chargers. When selecting the best candidates, a charging station utilisation factor should also be considered (Figure 7b). The final end stations targeted for electrification are then obtained by taking the cross-section of these two criteria (end stations underlined in red in Figure 7).     Figure 5 confirms that the share of total time of buses being parked in the depot is relatively small (approx. 30%) and comparable to the share of end-station parking time (approx. 25%). The rest of the time the buses spend in driving (40%), while only for a small portion of time (approx. 5%), they rest elsewhere, typically at bus stops. Figure 7a indicates that there are significant differences in bus resting time at different end stations. The average resting durations for most pronounced end stations are between 10 and 20 min, thus making them good candidates for installation of fast chargers. When selecting the best candidates, a charging station utilisation factor should also be considered (Figure 7b). The final end stations targeted for electrification are then obtained by taking the cross-section of these two criteria (end stations underlined in red in Figure 7).  Figure 6 indicates that the average bus velocity when operating is around 30 km/h and it is higher in early morning and night hours, as well as over weekends. The low-velocity gap between 2 a.m. and 5 a.m. corresponds to the interval when most of the buses are parked in the depot. The buses rarely visit the depot in other time intervals, particularly over the work days (when the average number of buses is lower than 1). Figure 5 confirms that the share of total time of buses being parked in the depot is relatively small (approx. 30%) and comparable to the share of end-station parking time (approx. 25%). The rest of the time the buses spend in driving (40%), while only for a small portion of time (approx. 5%), they rest elsewhere, typically at bus stops. Figure 7a indicates that there are significant differences in bus resting time at different end stations. The average resting durations for most pronounced end stations are between 10 and 20 min, thus making them good candidates for installation of fast chargers. When selecting the best candidates, a charging station utilisation factor should also be considered (Figure 7b). The final end stations targeted

General Description
As illustrated in flowchart in Figure 8, the EBSM simulates different types of city buses (CONV, HEV, PHEV and BEV) over the recorded driving cycles extracted by the DPPM (Section 3). The simulation first involves loading of vehicle-related parameters from the database, which need to be previously defined in the DMM. Next, the vehicle is simulated over the selected driving cycles by using the numerically-efficient backward-looking model (Subsection 4.2). Note that the vehicle model includes a control strategy that manages the gear ratio in the CONV and BEV cases, and the internal combustion engine (ICE) torque in the HEV and PHEV cases (Subsection 4.3). The emphasis has been on transforming the previously developed control strategy [17,18] to a form of off-line optimised maps, instead of using an on-line optimisation algorithm. The EBSM outputs time responses of key powertrain variables, such as cumulative fuel and electricity consumption, CO2 emissions and transmission gear ratio.

General Description
As illustrated in flowchart in Figure 8, the EBSM simulates different types of city buses (CONV, HEV, PHEV and BEV) over the recorded driving cycles extracted by the DPPM (Section 3). The simulation first involves loading of vehicle-related parameters from the database, which need to be previously defined in the DMM. Next, the vehicle is simulated over the selected driving cycles by using the numerically-efficient backward-looking model (Section 4.2). Note that the vehicle model includes a control strategy that manages the gear ratio in the CONV and BEV cases, and the internal combustion engine (ICE) torque in the HEV and PHEV cases (Section 4.3). The emphasis has been on transforming the previously developed control strategy [17,18] to a form of off-line optimised maps, instead of using an on-line optimisation algorithm. The EBSM outputs time responses of key powertrain variables, such as cumulative fuel and electricity consumption, CO 2 emissions and transmission gear ratio.

General Description
As illustrated in flowchart in Figure 8, the EBSM simulates different types of city buses (CONV, HEV, PHEV and BEV) over the recorded driving cycles extracted by the DPPM (Section 3). The simulation first involves loading of vehicle-related parameters from the database, which need to be previously defined in the DMM. Next, the vehicle is simulated over the selected driving cycles by using the numerically-efficient backward-looking model (Subsection 4.2). Note that the vehicle model includes a control strategy that manages the gear ratio in the CONV and BEV cases, and the internal combustion engine (ICE) torque in the HEV and PHEV cases (Subsection 4.3). The emphasis has been on transforming the previously developed control strategy [17,18] to a form of off-line optimised maps, instead of using an on-line optimisation algorithm. The EBSM outputs time responses of key powertrain variables, such as cumulative fuel and electricity consumption, CO2 emissions and transmission gear ratio.

Considered City Buses
The MAN Lion's City buses with the length of 12 m and the capacity of up to 126 passengers are represented in virtual simulation by the Volvo 7900 bus model. The Volvo 7900 platform was chosen because it includes all three e-bus variants considered (HEV, PHEV and BEV; Table 2). Volvo e-buses use lithium iron phosphate (LFP) battery due to its high specific power required for propulsion and fast charging. The battery packs of Volvo 7900 HEV, PHEV and BEV bus variants have energy capacities of 4.8 kWh, 19 kWh and 76 kWh, respectively ( Table 2).

Modelling
In the backward-looking models, the powertrain variables are calculated in the direction from the wheels towards the engine and/or e-motor, starting from the wheel speed and torque being defined by the driving cycles [20]. In order to boost the computational efficiency, the powertrain dynamics is neglected, except for the battery state-of-charge (SoC) dynamics that are represented by a first-order model.
The considered parallel configuration of a HEV/PHEV-type bus is illustrated in Figure 9a [17]. The battery is represented by the equivalent battery circuit model shown in Figure 9b, which is described by the following state equation [20,21]: where U oc is the open-circuit voltage, R is the internal resistance, I batt is the battery current, Q max is the maximal battery charge capacity and the SoC is defined as SoC = Q/Q max , with Q denoting the actual charge.

Considered City Buses
The MAN Lion's City buses with the length of 12 m and the capacity of up to 126 passengers are represented in virtual simulation by the Volvo 7900 bus model. The Volvo 7900 platform was chosen because it includes all three e-bus variants considered (HEV, PHEV and BEV; Table 2).
Volvo e-buses use lithium iron phosphate (LFP) battery due to its high specific power required for propulsion and fast charging. The battery packs of Volvo 7900 HEV, PHEV and BEV bus variants have energy capacities of 4.8 kWh, 19 kWh and 76 kWh, respectively ( Table 2). In the backward-looking models, the powertrain variables are calculated in the direction from the wheels towards the engine and/or e-motor, starting from the wheel speed and torque being defined by the driving cycles [20]. In order to boost the computational efficiency, the powertrain dynamics is neglected, except for the battery state-of-charge (SoC) dynamics that are represented by a first-order model.
The considered parallel configuration of a HEV/PHEV-type bus is illustrated in Figure 9a [17]. The battery is represented by the equivalent battery circuit model shown in Figure 9b, which is described by the following state equation [20,21]: where Uoc is the open-circuit voltage, R is the internal resistance, Ibatt is the battery current, Qmax is the maximal battery charge capacity and the SoC is defined as SoC = Q / Qmax, with Q denoting the actual charge. The model input P batt represents the battery output power defined as: where τ MG and ω MG are the motor/generator (M/G) torque and speed, respectively, η MG is the M/G machine efficiency (represented by η MG (ω MG ,τ MG ) map; see Figure 10c and [17]) and the coefficient k is equal to 1 or −1 depending on whether the M/G machine operates as a generator or motor, respectively. The M/G machine speed and torque are given by the following kinematic equations: where v v is the vehicle velocity, r w is the tire effective radius, h and i 0 are the transmission and final drive ratios, respectively, τ w and w w are the total wheels torque and speed, respectively, η tr (τ w ) and P 0 (ω w ) are drivetrain efficiency and idle power loss maps [17] and τ e is the engine torque considered as a control variable (in addition to h). The wheel torque is determined according to vehicle longitudinal dynamics equation covering the vehicle acceleration torque and aerodynamic, road grade and rolling resistances [20].
Energies 2020, 13, x FOR PEER REVIEW 9 of 25 The model input Pbatt represents the battery output power defined as: where τMG and ωMG are the motor/generator (M/G) torque and speed, respectively, ηMG is the M/G machine efficiency (represented by ηMG(ωMG,τMG) map; see Figure 10c and [17]) and the coefficient k is equal to 1 or −1 depending on whether the M/G machine operates as a generator or motor, respectively. The M/G machine speed and torque are given by the following kinematic equations: where vv is the vehicle velocity, rw is the tire effective radius, h and i0 are the transmission and final drive ratios, respectively, τw and ww are the total wheels torque and speed, respectively, ηtr(τw) and P0(ωw) are drivetrain efficiency and idle power loss maps [17] and τe is the engine torque considered as a control variable (in addition to h). The wheel torque is determined according to vehicle longitudinal dynamics equation covering the vehicle acceleration torque and aerodynamic, road grade and rolling resistances [20]. The fuel consumption at the driving cycle end time tf is determined as: The fuel consumption at the driving cycle end time t f is determined as: where A ek is the engine specific fuel consumption given by the map shown in Figure 11, ρ fuel is the diesel fuel density (ρ fuel = 845 g/L) and ω e equals ω MG or 0 when the engine is switched on or off, respectively. Note that the integral Equation (5) is realised by using the Euler integration method with the common sample time of the backward model equal to 1 s. where Aek is the engine specific fuel consumption given by the map shown in Figure 11, ρfuel is the diesel fuel density (ρfuel = 845 g/L) and ωe equals ωMG or 0 when the engine is switched on or off, respectively. Note that the integral Equation (5) is realised by using the Euler integration method with the common sample time of the backward model equal to 1 s. In the case of a conventional (CONV) bus, the battery and M/G machine are omitted from the functional scheme in Figure 9, while the AMT is replaced by a torque converter AT. The torque converter is represented by a backward-looking map ωι (ωt,τt) derived offline from the well-known static torque converter model (see [22] and references therein) and the map the subscripts i and t denote impeller/engine and turbine/transmission input variables, and Rτ (.) is the static model torque ratio map. For the BEV-type bus, the engine is omitted and a two-speed AMT is used.
Vehicle auxiliary devices (HVAC system, servo steering, air compressor, engine cooling fan and alternator) are modelled based on the nominal power of each device and a binary power-modulating signal, whose duty cycle is made dependent on the driving and atmospheric conditions (urban driving conditions and ambient temperature dependence are assumed) [23].

Control Strategy
In the BEV case, the aim of control strategy is to find an optimal gear ratio h * in each sampling instant, which minimises the battery discharging power and maximises the battery recharging/regenerative power: where ηbatt is the battery efficiency and Pd is the transmission input power demand calculated from vv and τw as shown in Figure 12. The optimal gear ratio h * is calculated offline and mapped as h * (Pd, ωw).
The SoC dependence of Pbatt has a minor effect on h * and is neglected in Eq. (6), and further on. Similarly, in the CONV case, the control strategy finds an optimal gear ratio h * * that minimises the fuel mass flow ṁf: The off-line obtained optimal solutions are mapped as h * (Pd, ωw).
The HEV and PHEV control strategy determines the transmission gear ratio h and the engine torque τe based on combining a rule-based (RB) controller and an equivalent consumption minimisation strategy (ECMS) [17,24]. The original RB+ECMS strategy is modified here to allow for computationally-efficient map-based realization, which is shown in Figure 12. A proportional-like In the case of a conventional (CONV) bus, the battery and M/G machine are omitted from the functional scheme in Figure 9, while the AMT is replaced by a torque converter AT. The torque converter is represented by a backward-looking map ω ι (ω t ,τ t ) derived offline from the well-known static torque converter model (see [22] and references therein) and the map τ i = τ t / R τ (ω t /ω i ), where the subscripts i and t denote impeller/engine and turbine/transmission input variables, and R τ (.) is the static model torque ratio map. For the BEV-type bus, the engine is omitted and a two-speed AMT is used.
Vehicle auxiliary devices (HVAC system, servo steering, air compressor, engine cooling fan and alternator) are modelled based on the nominal power of each device and a binary power-modulating signal, whose duty cycle is made dependent on the driving and atmospheric conditions (urban driving conditions and ambient temperature dependence are assumed) [23].

Control Strategy
In the BEV case, the aim of control strategy is to find an optimal gear ratio h * in each sampling instant, which minimises the battery discharging power and maximises the battery recharging/ regenerative power: where η batt is the battery efficiency and P d is the transmission input power demand calculated from v v and τ w as shown in Figure 12. The optimal gear ratio h* is calculated offline and mapped as h * (P d , ω w ). The SoC dependence of P batt has a minor effect on h* and is neglected in Equation (6), and further on. Similarly, in the CONV case, the control strategy finds an optimal gear ratio h* that minimises the fuel mass flowṁ f : The off-line obtained optimal solutions are mapped as h * (P d , ω w ).
Energies 2020, 13, 3410 11 of 23 dependent weighting factor W(eSoC) [17,24]: where W(.) is an SoC control error-dependent weighting function illustrated in Figure 12. For small SoC control errors, the 2D-ECMS is dominant, while for large SoC control errors, the RB+1D-ECMS is preferred to satisfy the engine power demand P * e (and, thus, the battery power demand P * batt). Finally, gear shift delay logic is implemented to prevent frequent gear switching [17]. The goal is to prevent gear shift occurrence in k th time step, i.e., rather use the gear ratio from the previous (k-1) th step, hk-1, if the time elapsed since the last gear shift tsh is lower than an arbitrarily set shift delay The HEV and PHEV control strategy determines the transmission gear ratio h and the engine torque τ e based on combining a rule-based (RB) controller and an equivalent consumption minimisation strategy (ECMS) [17,24]. The original RB+ECMS strategy is modified here to allow for computationally-efficient map-based realization, which is shown in Figure 12. A proportional-like battery SoC controller commands the battery recharging power-P batt , which is added to the transmission input power demand P d to obtain the engine power demand P * e . The demanded engine power P * e is compared with engine on and off thresholds P on and P off < P on , respectively, in order to determine the engine on/off state EN st . The engine is exceptionally kept switched on in the case P * e < P off when the speed-dependent M/G machine power limit is not high enough to satisfy the driver power demand P d in the fully electric driving mode. If the engine is switched on, the signal P * e is fed to the ECMS to find the optimal values of h and τ e . Otherwise, the electric driving mode is activated and the gear ratio is determined according to Equation (6).
In the original RB+ECMS strategy [17], the equivalent fuel consumptionṁ eq (P e , P d , ω w , h, SoC) is minimised instantaneously and on-line with respect to both control variables h and τ e . In the simplified map-based RB+ECMS version considered here, the equivalent fuel consumption is minimised in two stages. In the first stage, the ECMS is applied to discrete operating points along the constant power curve P * e (denoted in Figure 11 by blue circles) to determine the optimal gear ratio: whereṁ eq is the equivalent fuel consumption rate containing the actual fuel consumption rateṁ f and a battery power-equivalent fuel rate (see [17,24] for details). The off-line obtained optimal solutions are stored in a three-dimensional (3D) map h * (P * e , P d , ω w ) representing the RB+1D-ECMS control map.
In the second stage, the ECMS is applied along the engine torque axis (see green arrows in Figure 11): The off-line optimisation results are stored in a 2D map τ * e,2D (P d , ω e ) representing the 2D-ECMS control map. Finally, the engine torque obtained by the RB+1D-ECMS (as τ * e,RB = P * e /(h * i o ω w )) and the one obtained by the 2D-ECMS are combined/blended on-line using the SoC control error (e SoC )-dependent weighting factor W(e SoC ) [17,24]: where W(.) is an SoC control error-dependent weighting function illustrated in Figure 12. For small SoC control errors, the 2D-ECMS is dominant, while for large SoC control errors, the RB+1D-ECMS is preferred to satisfy the engine power demand P * e (and, thus, the battery power demand P * batt ). Finally, gear shift delay logic is implemented to prevent frequent gear switching [17]. The goal is to prevent gear shift occurrence in k th time step, i.e., rather use the gear ratio from the previous (k−1) th step, h k−1 , if the time elapsed since the last gear shift t sh is lower than an arbitrarily set shift delay threshold t th and if h k−1 gives feasible set u k of engine and M/G machine operating points in the k th step (denoted by Π): The same gear shifting delay logic is applied in the CONV and BEV cases. The shift delay threshold is set here to t th = 2 s.
The above described simplified, map-based control strategy has been found to result in a negligible model response deviation when compared to the use of original strategy, which was proven to be close to the dynamic programming-based global optimum [17]. On the other hand, the execution time is reduced by around 200 times. The achieved execution time, expressed as the amount of microseconds needed to simulate one second of real time (for a workstation having 16 GB RAM and Intel ® Xeon ® Processor E5-1620 v3 @ 3.50GHz) falls in the range from 50 to 87 µs/s depending on vehicle type. This results in approximate yearly 10-bus fleet simulation time ranging from approximately 4.5 h to 7.5 h, which is deemed acceptable for such a large-scale fleet simulation. Note that the execution time could further be reduced by using parallel computing.
The PHEV can operate in two characteristic modes [20]: (i) charge depleting (CD) followed by charge sustaining (CS), where the former involves the engine only when absolutely needed and the latter correspond to hybrid operation at the target SoC of 30%; (ii) blended mode where engine is regularly used all over the driving cycle for additional energy savings. For the sake of simplicity, the CD/CS mode is considered in this paper.

Simulation Results
The results related to relative fuel and/or electricity consumptions for different city bus types are given in Table 3 for the full recording period. The relative difference between the simulated (Sim) and recorded (Rec) fuel consumptions for the CONV bus is equal to only 1.4%. Therefore, the CONV simulation model used as a basis for e-bus modelling can be considered accurate. Note that although the real and simulated buses are different (MAN Lion City and Volvo 7900), the validation is considered fair, as the two buses are rather comparable in terms of size, mass, engine power, number of passengers and other similar factors.
The simulated electricity consumptions of PHEV-and BEV-type buses are close to recorded ones documented in the ZeEUS project report [25] for Volvo 7900 bus series (Table 3). In the PHEV case, the simulated fuel consumption is by 30% higher than the ZeEUS recorded one, but this discrepancy is compensated for by 26% higher recorded electricity consumption when compared to the simulated one. In the BEV case, the relative difference in electricity consumption equals 6%. The simulated HEV fuel consumption is reduced by 50% when compared to CONV simulation results, while the manufacturer states the fuel consumption reduction from 39% to 45% reported by operators [19]. The observed, relatively modest discrepancies in fuel/electricity consumption may be related to difference in considered bus weights (passenger weight is fixed to 1250 kg), road slope and traffic congestion conditions, as well as regenerative braking capacity (set to the maximum amount of 100% in simulation). * Recorded PHEV and BEV fuel and electricity consumptions are taken from the ZeEUS project report [25]. ** Estimated based on information on fuel consumption reduction from 39% to 45% for Volvo 7900 Hybrid vs. Volvo 7900 according to [19]. Figure 10 shows the PHEV-case simulation results for the recorded driving cycle shown in Figure 4 repeated 15 times and the initial battery SoC equal to 90%. When the CD mode is active, the engine is used only when needed and the cumulative fuel consumption is often constant (a stepwise-like response, Figure 10a). After entering the CS mode, the engine is more active to sustain the battery SoC (Figure 10b). The control strategy deploys the operating points of engine (when switched on) and M/G machine in the high efficiency areas of corresponding maps (Figure 10c-d), thus minimising the energy consumption.

General Description
A generic framework for virtual simulation of an e-bus fleet over the recorded driving cycles is represented by the flowchart shown in Figure 13. The model of city bus of any type (Section 4) is initialised based on the data stored in DMM and virtually run over the recorded driving cycles for the specified period of time, thus resulting in fuel and/or electricity consumption output data. In the case of a PHEV-or BEV-type bus the user needs to specify locations and types of charging stations, the nominal vehicle battery capacity and grid power constraints to execute the simulation. The simulation is repeatedly run for a peak day or peak week for different charging infrastructure and battery capacity specifications, in order to find a nearly-optimal configuration, which would be finally re-run for a given, longer period of time to calculate the fuel/electricity consumption and charging station utilisation statistics. In the case of BEV-type bus, the COM automatically adds reserve bus(es) if needed, and calculate their final number and related statistics. The bus intervening algorithm monitors each BEV's battery SoC, and when it drops below a predefined minimum value (0.3, herein), the immediate replacement with reserve e-bus is conducted. At the same time, it is taken into account that the depleted bus needs some constant time to reach the nearest charging station, where it is to be charged (as any other bus), and once it is fully charged, it will be waiting for the next replacement/intervention.

Charging Management Algorithm
Charging management is described by the flowchart shown in Figure 14. First, it is checked if a PHEV-or BEV-type bus has arrived to an end station/depot and if that station has a charger installed. If the charger is not occupied or if the bus has a lower battery SoC of any of the buses already being charged, the bus is put on charge; otherwise, it remains in the charging queue. Note that each station can be set to have an arbitrary number of chargers, as described with Figure 13.
reserve bus(es) if needed, and calculate their final number and related statistics. The bus intervening algorithm monitors each BEV's battery SoC, and when it drops below a predefined minimum value (0.3, herein), the immediate replacement with reserve e-bus is conducted. At the same time, it is taken into account that the depleted bus needs some constant time to reach the nearest charging station, where it is to be charged (as any other bus), and once it is fully charged, it will be waiting for the next replacement/intervention.

Charging Management Algorithm
Charging management is described by the flowchart shown in Figure 14. First, it is checked if a PHEV-or BEV-type bus has arrived to an end station/depot and if that station has a charger installed. If the charger is not occupied or if the bus has a lower battery SoC of any of the buses already being sampling instant k according to: where SoCk is the current SoC and SoCf is the target SoC at departure. Inserting Equation (12) into the battery state Equation (1) yields the charging power Pbatt < 0 to be applied in the k th sampling instant: If the charging power −Pbatt calculated from Equation (13) is greater/less than the maximum/minimum allowable power (defined by the charger selected), the charging power is limited to the maximum/minimum power, respectively. Note that ∆ = − is saturated in Equation (13) to its lower limit of 30 s to avoid division by zero.
Once the charging power profile is obtained for each sampling step k and for i th vehicle from the total number of Nv vehicles connected to chargers at the same grid sections (e.g., depot), it is checked if the total charging power is greater than the maximum grid power , . If this applies, the charging power is scaled down to satisfy the grid power constraint:

Obtaining of Near-Optimal Charging System Configurations
According to the city bus transport characterisation results from Figure 7, there is a number of end stations with relatively long bus resting durations and potentially high utilization of charger units. Additionally, the end station resting time share approaches that of depot ( Figure 5), and there are no other emphasised stop locations. Therefore, fast charging stations and belonging transformer substations can be installed at end stations to provide bus recharging, while otherwise the available The charging process is managed by taking into account the requirements on satisfying the departure schedule, minimising battery power loss and respecting the grid power constraints. According to [26,27], the battery energy loss is minimised by demanding a linear change of SoC all over the remaining charging interval ∆T ch = t f − t k . Therefore, the SoC rate is updated in each sampling instant k according to: where SoC k is the current SoC and SoC f is the target SoC at departure. Inserting Equation (12) into the battery state Equation (1) yields the charging power P batt < 0 to be applied in the k th sampling instant: If the charging power −P batt calculated from Equation (13) is greater/less than the maximum/minimum allowable power (defined by the charger selected), the charging power is limited to the maximum/minimum power, respectively. Note that ∆T ch = t f − t k is saturated in Equation (13) to its lower limit of 30 s to avoid division by zero. Once the charging power profile is obtained for each sampling step k and for i th vehicle from the total number of N v vehicles connected to chargers at the same grid sections (e.g., depot), it is checked if the total charging power is greater than the maximum grid power P grid,max . If this applies, the charging power is scaled down to satisfy the grid power constraint:

Obtaining of Near-Optimal Charging System Configurations
According to the city bus transport characterisation results from Figure 7, there is a number of end stations with relatively long bus resting durations and potentially high utilization of charger units. Additionally, the end station resting time share approaches that of depot ( Figure 5), and there are no other emphasised stop locations. Therefore, fast charging stations and belonging transformer substations can be installed at end stations to provide bus recharging, while otherwise the available power can be utilised to supply city e-mobility hubs built around the end stations. High-power off-board chargers with built-in pantograph are considered (150 or 300 kW, see Table 2) [19]. In addition, the slow-to-modestly fast plug-in charging solutions can be considered for a depot, where the charging time can be long in night ( Figure 6). Figure 15 shows the PHEV fleet simulation results for different number of end stations equipped with a single fast charger per station (150 kW) and a five work day period. Charging in depot was not considered because the small-capacity PHEV battery (Table 2) can quickly be recharged at the end stations, where the buses rest for a relatively long time (Figures 6 and 7). The results shown in Figure 15 point out that the fuel consumption saving converges to −41% as the number of end station charging spots approaches six. Of course, as the fuel consumption reduces, the electricity consumption grows, but the overall energy cost is reduced by 17% due to cheaper electricity. By conducting PHEV fleet simulations over the five-month period, it has been found that the optimal number of charging stations should be incremented to seven.  Table 2) [19]. In addition, the slow-to-modestly fast plug-in charging solutions can be considered for a depot, where the charging time can be long in night ( Figure 6). Figure 15 shows the PHEV fleet simulation results for different number of end stations equipped with a single fast charger per station (150 kW) and a five work day period. Charging in depot was not considered because the small-capacity PHEV battery (Table 2) can quickly be recharged at the end stations, where the buses rest for a relatively long time (Figures 6 and 7). The results shown in Figure 15 point out that the fuel consumption saving converges to -41% as the number of end station charging spots approaches six. Of course, as the fuel consumption reduces, the electricity consumption grows, but the overall energy cost is reduced by 17% due to cheaper electricity. By conducting PHEV fleet simulations over the five-month period, it has been found that the optimal number of charging stations should be incremented to seven.

Case of BEV Fleet
The BEV fleet simulation results are shown in Table 4. The full five-month period is considered to cover a larger number of "critical" days when reserve buses may be needed. Only scenarios with the number of end-station charging spots being in the vicinity of the optimal one found for the PHEV fleet is examined (around six stations plus depot, each considering a single fast charger). The maximum charging power is set to the levels of 150 kW or 300 kW. Finally, various battery capacities are considered (76, 150 and 250 kWh, as the capacities available for the considered bus [19]).

Case of BEV Fleet
The BEV fleet simulation results are shown in Table 4. The full five-month period is considered to cover a larger number of "critical" days when reserve buses may be needed. Only scenarios with the number of end-station charging spots being in the vicinity of the optimal one found for the PHEV fleet is examined (around six stations plus depot, each considering a single fast charger). The maximum charging power is set to the levels of 150 kW or 300 kW. Finally, various battery capacities are considered (76, 150 and 250 kWh, as the capacities available for the considered bus [19]). Table 4. BEV fleet simulation results for different number of charging spots (located at end stations and depot) and reserve buses, and different battery capacities (full five-month period). The results shown in Table 4 point out that by increasing the number of charging stations, the percentage of total electricity consumed by reserve buses drops from 9.2% (case BEV 1) to 1.8% (case BEV 2). Likewise, the number of bus swaps (concerning reserve bus) drops from 558 in 106 (out of 152) days (BEV 1) to 94 in 54 days (BEV 2). Figure 16 indicates that in the case BEV 2 notable bus swaps occurs only in several days, which are characterised by peak traffic load (typically due to specific needs such as moving tourists from cruising ships to the old city). Similar trends apply to the case of increasing the battery capacity from 76 kWh to 150 kWh and further to 250 kWh (cases BEV 4 and BEV 5, respectively), where the reserve buses are marginally needed in the former case, and not needed in the latter case. The results shown in Table 4 point out that by increasing the number of charging stations, the percentage of total electricity consumed by reserve buses drops from 9.2% (case BEV 1) to 1.8% (case BEV 2). Likewise, the number of bus swaps (concerning reserve bus) drops from 558 in 106 (out of 152) days (BEV 1) to 94 in 54 days (BEV 2). Figure 16 indicates that in the case BEV 2 notable bus swaps occurs only in several days, which are characterised by peak traffic load (typically due to specific needs such as moving tourists from cruising ships to the old city). Similar trends apply to the case of increasing the battery capacity from 76 kWh to 150 kWh and further to 250 kWh (cases BEV 4 and BEV 5, respectively), where the reserve buses are marginally needed in the former case, and not needed in the latter case.

Case
Based on the above results, the case BEV 5 might be considered as optimal. However, since increasing of the battery capacity of each bus in a fleet is rather expensive and the need for reserve buses in case BEV 2 is minor (only 1.8%), the case BEV 2 has been adopted as an optimal for final simulations discussed in Subsection 5.4.  Table 5 gives the energy consumption simulation results for the PHEV and BEV fleets configured in the previous two sections, as well as HEV and CONV fleets, all simulated over the full recording period. The HEV fleet can reduce the fuel consumption by around 50% compared to CONV fleet, owing to regenerative braking, switching the engine off in low speed conditions and placing the engine operating points in the high-efficiency region (Section 4; see also Table 3). Using the PHEV buses provides additional fuel savings, i.e., the fuel consumption reduction compared to CONV case is around 70%. This is due to the use of electricity coming from grid by means of fast charging at end stations. The PHEV fleet electricity consumption equals almost 50% of what is consumed by BEV fleet, which is owing to a relatively low portion of operation in the CS mode (25.5%).

Comparative Energy Consumption Results
The COM also calculates the well-to-wheel CO2 emissions taking into account the data on CO2 emissions of power plants, as used in [28] for the cases considering coal, natural gas and renewable energy production. The HEV vs. CONV fleet CO2 emissions reduction is around 50%, while in the cases of PHEV and BEV fleets, the reduction is from 30% to 65% and from 30% to 93%, respectively, where the lower and higher margins correspond to coal and renewable energy production scenarios. Table 5. Fuel and electricity consumption per fleet type for full five-month period. Based on the above results, the case BEV 5 might be considered as optimal. However, since increasing of the battery capacity of each bus in a fleet is rather expensive and the need for reserve buses in case BEV 2 is minor (only 1.8%), the case BEV 2 has been adopted as an optimal for final simulations discussed in Section 5.4. Table 5 gives the energy consumption simulation results for the PHEV and BEV fleets configured in the previous two sections, as well as HEV and CONV fleets, all simulated over the full recording period. The HEV fleet can reduce the fuel consumption by around 50% compared to CONV fleet, owing to regenerative braking, switching the engine off in low speed conditions and placing the engine operating points in the high-efficiency region (Section 4; see also Table 3). Using the PHEV buses provides additional fuel savings, i.e., the fuel consumption reduction compared to CONV case is around 70%. This is due to the use of electricity coming from grid by means of fast charging at end stations. The PHEV fleet electricity consumption equals almost 50% of what is consumed by BEV fleet, which is owing to a relatively low portion of operation in the CS mode (25.5%). The COM also calculates the well-to-wheel CO 2 emissions taking into account the data on CO 2 emissions of power plants, as used in [28] for the cases considering coal, natural gas and renewable energy production. The HEV vs. CONV fleet CO 2 emissions reduction is around 50%, while in the cases of PHEV and BEV fleets, the reduction is from 30% to 65% and from 30% to 93%, respectively, where the lower and higher margins correspond to coal and renewable energy production scenarios.

General Description and TCO Model
According to Figure 2, the TEAM uses the simulation data outputted by the COM, as well as the fleet loan payment, insurance, registration, maintenance and similar costs provided by DMM, in order to calculate the fleet TCO ( Figure 17). The TCO corresponds to what is in financial terminology called Net Present Value (NPV) of an investment, which is an index that valorises the investment while considering the time value of money. Rates at which the money value decreases or increases over time are in this case modelled by the inflation and discount rates, respectively (where the latter corresponds to the profit that today's money can generate in the future through investments or bank savings). Calculation of future value of money is called compounding, while the opposite approach, in which the NPV of future money is calculated, is referred as discounting. The TCO is calculated by discounting all future expenses, which the investment is expected to generate, to the present time, as shown in Figure 17.
The TCO model components ( Figure 17) are divided into three groups depending on the time basis on which the input expenses data are sampled [29], and the corresponding individual costs are given in Table 6 (with no VAT included). The bus service life is considered to be 12 years, the inflation rate 3% and the discount rate 7%. The annually sampled data include registration, maintenance and insurance (RMI) cost, which have been determined for the CONV fleet based on the (past) data provided by the city bus transport operator, and discounted to prices in 2019 according to inflation data [30]. The RMI cost for the HEV, PHEV and BEV fleets are assumed to be 15%, 20% and 40% lower, respectively, when compared to the CONV fleet, because of the significantly reduced CO 2 emissions and simplified maintenance of e-buses [31][32][33].
The monthly expenses relate to loan payment for purchase of new vehicles and charging infrastructure, including the cost of replacing the e-bus batteries. A general-purpose bank loan is assumed, which is taken over a period of seven years, with a continuous interest rate of 5% and equal monthly annuities. The daily sampled data relate to operating cost, i.e., the fuel and electricity expenses, which are calculated by multiplying fuel and/or electricity consumptions obtained by COM simulations with fuel and/or electricity prices. As in the case of the annually sampled data, the operational cost is adjusted for inflation. Irregular maintenance cost is modelled by a fixed rate occurring every two years. The TEAM also provides the possibility of sensitivity analysis, which allows for the investigation of to what extent variations of a particular parameter affect the TCO. This helps to determine the TCO model reliability, e.g., parameters that cause higher TCO sensitivity should be more reliably estimated. The sensitivity analysis is not considered in this paper.
fleet loan payment, insurance, registration, maintenance and similar costs provided by DMM, in order to calculate the fleet TCO ( Figure 17). The TCO corresponds to what is in financial terminology called Net Present Value (NPV) of an investment, which is an index that valorises the investment while considering the time value of money. Rates at which the money value decreases or increases over time are in this case modelled by the inflation and discount rates, respectively (where the latter corresponds to the profit that today's money can generate in the future through investments or bank savings). Calculation of future value of money is called compounding, while the opposite approach, in which the NPV of future money is calculated, is referred as discounting. The TCO is calculated by discounting all future expenses, which the investment is expected to generate, to the present time, as shown in Figure 17. The TCO model components ( Figure 17) are divided into three groups depending on the time basis on which the input expenses data are sampled [29], and the corresponding individual costs are given in Table 6 (with no VAT included). The bus service life is considered to be 12 years, the inflation rate 3% and the discount rate 7%. The annually sampled data include registration, maintenance and insurance (RMI) cost, which have been determined for the CONV fleet based on the (past) data provided by the city bus transport operator, and discounted to prices in 2019 according to inflation data [30]. The RMI cost for the HEV, PHEV and BEV fleets are assumed to be 15%, 20% and 40% lower, respectively, when compared to the CONV fleet, because of the significantly reduced CO2 emissions and simplified maintenance of e-buses [31][32][33].
The monthly expenses relate to loan payment for purchase of new vehicles and charging infrastructure, including the cost of replacing the e-bus batteries. A general-purpose bank loan is assumed, which is taken over a period of seven years, with a continuous interest rate of 5% and equal   [34]. 2 Costs for transformer substation (TS) and charging station (CS) are estimated based on the data provided by local electric utility company and [35], respectively. 3 The battery replacement costs are estimated based on [35] and the replacement is assumed to occur every 6 years because the average bus battery life is 5−12 years [36]. 4

Simulation Results
The TCO results are given in Figure 18 for different types of city bus fleets and charging configurations selected in Section 5. Different charging scenarios are considered, starting from optimistic Scenario 1 to conservative Scenario 6. In the basic case (Scenario 1), all the e-bus fleets turn out to be competitive with the CONV fleet, which is explained by the influence of high share of fuel cost (see Figure 19) for the particular case of relatively large fleet utilisation (250 km/bus/day in average; see also Figure 5). For the same reason, the e-bus fleets have relatively comparable TCO values. Similar results are obtained if the fuel and electricity prices are randomly sampled (Scenario 2), rather than being constant as given in Table 6. Scenario 3 accounts for the need to use reserve buses in the case of BEV fleet, as found by COM simulations (Section 5). Since in the considered case, BEV 2, the use of two reserve buses results in a marginal increase of electricity consumption (Table 4) and a low number of bus swaps, the use of second reserve bus is very marginal, and is thus excluded from the TCO analysis. Due to the cost of reserve bus, the BEV fleet TCO increases above that of PHEV fleet, but it is still competitive to CONV fleet. When accounting for the e-buses' battery replacement cost (Scenario 4), the TCO of BEV fleet, which has the largest and costliest battery, becomes around 10% higher than that of CONV fleet. If the PHEV-and BEV-type bus electricity consumption is increased by the factor of 40% (Scenario 5) or 100% (Scenario 6) to account for modelling errors (e.g., those related to heating system in winter), the PHEV fleet becomes marginally competitive or uncompetitive, respectively, while the BEV vs. CONV fleet TCO excess tops 23%. This TCO excess in the ultimate BEV case may be compensated for by larger incentives, higher ticket prices (which would reflect better passenger experience), future increase in fuel prices, future decrease of battery prices and similar factors. The comparative TCO time profiles for different types of bus fleet are shown in Figure 20 for Scenario 4, which is deemed to be most realistic scenario involving the battery replacement and reserve bus cost. The corresponding time profiles of individual TCO costs are shown in Figure 21. For the PHEV, and particularly the BEV fleet, the TCO rapidly rises during the first 7 years due to loan expenses related to the purchase of these expensive buses and corresponding charging infrastructure ( Figure 20). Once the loan is paid off, the energy cost becomes dominant, where the efficiency of e-buses and low cost of electricity become beneficial and bring significant savings, as opposed to the CONV case, where the fuel expenses dominate ( Figure 21). Figure 19 shows the percentage shares of individual costs for different types of bus fleets. As the electrification evolves from HEV, via PHEV to BEV buses, the energy (fuel and electricity) cost share monotonically and significantly reduces, but the bus and charging infrastructure cost share increases with similar trends. The PHEV and particularly BEV fleets have lower RMI cost, but this saving is not large enough to compensate for the battery replacement cost.    Figure 20 for Scenario 4, which is deemed to be most realistic scenario involving the battery replacement and reserve bus cost. The corresponding time profiles of individual TCO costs are shown in Figure 21. For the PHEV, and particularly the BEV fleet, the TCO rapidly rises during the first 7 years due to loan expenses related to the purchase of these expensive buses and corresponding charging infrastructure ( Figure 20). Once the loan is paid off, the energy cost becomes dominant, where the efficiency of e-buses and low cost of electricity become beneficial and bring significant savings, as opposed to the CONV case, where the fuel expenses dominate (Figure 21).

Conclusions
A unique city bus transport electrification tool has been presented, which relies on virtual simulation of conventional and e-bus fleets over recorded driving cycles, including charging management. The tool allows for arbitrary setting of charging station locations, number and parameters of chargers per each station, grid constraints, e-bus battery capacity and similar factors. A special emphasis was on providing a computationally efficient, multiple-maps-based backward looking e-bus model including the vehicle control strategy and heuristic charging management algorithm. In support of electrification planning and decision making, the tool finally provides the total cost of ownership (TCO) of city bus fleets with conventional (CONV) buses and different types of e-buses (HEV, PHEV and BEV). It has been demonstrated through a pilot study related to city bus

Conclusions
A unique city bus transport electrification tool has been presented, which relies on virtual simulation of conventional and e-bus fleets over recorded driving cycles, including charging management. The tool allows for arbitrary setting of charging station locations, number and parameters of chargers per each station, grid constraints, e-bus battery capacity and similar factors. A special emphasis was on providing a computationally efficient, multiple-maps-based backward looking e-bus model including the vehicle control strategy and heuristic charging management algorithm. In support of electrification planning and decision making, the tool finally provides the total cost of ownership (TCO) of city bus fleets with conventional (CONV) buses and different types of e-buses (HEV, PHEV and BEV). It has been demonstrated through a pilot study related to city bus transport in Dubrovnik, based on which the following main conclusions are drawn.
1) The considered city bus transport system is such that the city buses are resting in the depot  Figure 19 shows the percentage shares of individual costs for different types of bus fleets. As the electrification evolves from HEV, via PHEV to BEV buses, the energy (fuel and electricity) cost share monotonically and significantly reduces, but the bus and charging infrastructure cost share increases with similar trends. The PHEV and particularly BEV fleets have lower RMI cost, but this saving is not large enough to compensate for the battery replacement cost.

Conclusions
A unique city bus transport electrification tool has been presented, which relies on virtual simulation of conventional and e-bus fleets over recorded driving cycles, including charging management. The tool allows for arbitrary setting of charging station locations, number and parameters of chargers per each station, grid constraints, e-bus battery capacity and similar factors. A special emphasis was on providing a computationally efficient, multiple-maps-based backward looking e-bus model including the vehicle control strategy and heuristic charging management algorithm. In support of electrification planning and decision making, the tool finally provides the total cost of ownership (TCO) of city bus fleets with conventional (CONV) buses and different types of e-buses (HEV, PHEV and BEV). It has been demonstrated through a pilot study related to city bus transport in Dubrovnik, based on which the following main conclusions are drawn.
(1) The considered city bus transport system is such that the city buses are resting in the depot during a relatively short period over the night (typically 3 h), while they are dwelling at end stations for rather significant time (from 10 to 20 min per stay). Therefore, fast charging at end stations (and also in depot for BEV-type buses) relying on stationary chargers equipped with pantograph has been found to be a favourable solution. (2) The use of a specific, map-based structure of the e-bus model allowed for simulating the bus fleets 20,000 times faster than real time, thus, reducing the full-year 10-bus fleet simulation to a couple of hours on a standard computer workstation. (3) The comparative virtual simulation results have shown that the use of HEV-and PHEV-type city buses results in reduction of fuel consumption of up to 50% and 70%, respectively, when compared to CONV buses, while BEV buses do not consume fuel at all. The charging system optimisation has shown that the optimal number of end stations equipped with fast chargers is seven (out of 10), where a single reserve bus is marginally needed in the BEV case. The BEV battery capacity can be relatively small (76 kWh) due to the effective opportunity charging and relatively short routes. (4) The TCO analysis has pointed out that the BEV fleet cannot be competitive to CONV fleet (8.6% higher TCO for BEV vs. CONV), while the HEV fleet is competitive (12.8% lower TCO vs. CONV) and the PHEV fleet is marginally competitive (3.8% lower TCO vs. CONV) in a realistic scenario involving the battery replacement and single reserve bus in the BEV case (Scenario 4). Although the HEV fleet is competitive to the CONV fleet and can reduce the fuel consumption and emissions by up to 50%, it still shares the basic disadvantages of CONV fleet (noisy, no e-drive option in low emission zones, significant emissions).
The future work could be directed to the following tool improvements: (i) optimisation of charging system configuration (e.g., by using genetic algorithm) instead of using the expert knowledge when repeating the virtual simulations for different configurations; (ii) using parallel computing to further increase the numerical efficiency of virtual simulations, particularly when optimisation loop is added; (iii) off-line route rescheduling to avoid or mitigate the use of reserve buses in the BEV case or to redirect the e-buses to designated location for daily recharging; (iv) performing city bus transport electrification analyses using limited set of standard GPS/GPRS tracking data while relying more on the bus schedules; (v) considering other economic models when conducting techno-economic analyses to determine the TCO. Funding: It is gratefully acknowledged that this work was supported by the EU European Regional Development Fund trough Interreg CE project SOLEZ ("Smart Solutions supporting Low Emission Zones and other low-carbon mobility policies in EU cities"; https://www.interreg-central.eu/Content.Node/SOLEZ.html).
Acknowledgments: In addition, the research work of the first author has been supported by the European Regional Development Fund under the grant KK.01.1.1.01.0009 (DATACROSS). The authors' appreciation also goes to the city bus transport company Libertas Dubrovnik d.o.o. and the Dubrovnik Development Agency (DURA) for the data support and effective cooperation.

Conflicts of Interest:
The authors declare no conflict of interest.