Multi-Level Energy Management for Hybrid Electric Vehicles — Part I

The fuel economy of a hybrid electric vehicle (HEV) is improved, by taking the energy relevant system states into account in the energy management system (EMS). With an increasing number of states and decision variables, energy optimizing algorithms in the EMS can be prohibitive for real-time implementation. In part I of this work, a model-based, multi-level approach is taken to subdivide the original (large) optimization problem into computational efficient sub-problems, based on optimal control techniques using a preview. The resulting EMS solves the problem of power-split between engine and motor/generator, mode and gear switching including switching costs, with battery energy constraints. The superior energy efficiency of the multi-level EMS is simulated on a representative heavy duty drive cycle, where it saves 7.0% fuel, compared to a conventional vehicle, where the baseline EMS for the HEV saves 5.8%. In part II, real-world validation of the EMS is performed.


Introduction
Hybrid electric vehicles (HEVs) have emerged as a promising solution to reduce operational cost in commercial road transportation, while complying to increasingly stringent emission legislation.Since HEVs have more than one power converter, they offer additional control freedom, compared to conventional vehicles, which give opportunities for the energy management system (EMS) to decrease fuel consumption and emissions.The EMS has to consider the energy relevant systems of the HEV and already a large amount of solutions have been proposed that take the battery energy dynamics into account [1].However, the system efficiency can be further improved by taking additional system states into account as suggested by [2,3], referred to as a 'unified', 'integrated', 'total' or 'holistic' energy management.Examples of such additional systems are the battery with its temperature and aging characteristics [4,5], engine after-treatment system [3,6], waste-heat recovery system [7], combustion engine [8], or the cabin heater [9].
Using model information of these systems, the EMS can be posed as an optimization problem, in order to ensure energy efficiency.The computational complexity of this problem, however, increases with the number of states, presence of state constraints, number of decision variables, type of decision variables (continuous/discrete) and nonlinearity of the models.Many methods exist to design an optimizing EMS, using e.g., Dynamic Programming [10,11], Pontryagin's Minimum Principle [12] or Model Predictive Control [13,14] and generate close to optimal results; however they can be computationally too demanding for real-time implementation.In this work, a combination of methods is used, to ensure an efficient computation for each of the sub-problems.An overview of real-time implementable optimization approaches for a (parallel) hybrid vehicle with the decision variables: power-split, stop-start and gear selection, is given in Table 1.In particular, discrete control variables increase the computational complexity, explaining why gear selection is often omitted from the EMS or solved in a separate step with heuristics [15][16][17].In [18], an integrated approach is presented, where gear selection and stop-start are part of a sequential optimization with the power split.In these approaches, the cost of stop-start and gear change events are not considered, which can result in unacceptable switching behaviour, like hunting oscillations.To overcome this problem, costs on switchings are included, solved using Dynamic Programming (DP) and Quadratic Programming (QP) [19] or as one Mixed Integer Linear Program (MILP) [20].For real-time solving, the allowed model complexity is limiting and adding additional states to these EMSs is computationally prohibitive.By partitioning a large optimization problem into a set of smaller problems, as in a distributed control system [21], the computational efficiency and robustness are improved, albeit losing the guarantee of global optimality.The process of partitioning is not trivial and many different structures exist.For the power-split problem, two-level structures can be identified in literature [22][23][24], where planning of the battery energy is separated from the power-split decisions.In this work, a novel EMS is developed as a multi-level system [25,26], in which the higher layers have an increasing level of abstraction of the system, to be able to optimize and coordinate the lower layers [27], while having a large decoupling between the levels.The multi-level EMS solves this control problem with preview, for real-time implementation, using optimal control techniques (PMP and DP), thereby eliminating calibration of parameters.In Figure 1, the proposed control system is illustrated.Each level minimizes fuel, however, dependent on the function of the level, different model information and corresponding optimal control techniques are used.

•
On the first level, the power-split is explicitly solved using the Pontryagin Minimum Principle (PMP), starting from [15].This method is extended with costs on mode and gear switching, thereby eliminating unacceptable switching behaviour.A Dynamic Programming (DP) routine solves the discrete subproblem.Route and vehicle information determine the switching costs, supporting a full model based approach.

•
The second level optimizes the battery state of charge with input-and state-constraints, and provides efficiency information of the hybrid system to the other levels using the battery costate from the PMP solution.Due to the mode switching system on the vehicle, charge sustaining behaviour must be enforced, by using an additional switching algorithm between non-unique solutions.

•
The third level provides the necessary route information to the lower layers.By predicting the velocity along the route, using road slope and velocity limitations, the road load on the driveline is determined.In Section 2, the model of the parallel hybrid vehicle is described, after which the multi-level EMS for this vehicle is formulated in Section 3. The algorithms used on each of the three levels are explained in Section 4 (power-split, mode and gear switch optimization), Section 5 (battery energy optimization) and Section 6 (velocity prediction).Offline and online solution schemes of the multi-level optimization are compared in Section 7, where high fidelity simulation results show the fuel benefit of the algorithm.The conclusions are drawn in Section 8.The real-world validation of the multi-level EMS is described in [28].

Parallel Hybrid Electric Vehicle Model
For the design of a model-based EMS, the applied models of the hybrid electric vehicle (HEV) are described in this section.The parallel hybrid vehicle under consideration is schematically depicted in Figure 2. The fuel power P f flows from the tank (FT) to the Combustion Engine (CE), converted to mechanical power P e , thereby depleting the available fuel energy E f .Dependent on the position of the clutches CL 1 and CL 2 , the power from the Motor/Generator (MG) P m is added to P e , resulting in the power at the power-split point P p .This power is transferred through the gearbox, final drive and wheels (GB), resulting in the driveline power P d of the vehicle (VH).Dependent on the road load acting on the vehicle (F 0 ), the travelled distance s and velocity v will change.The MG exchanges electrical power P b with the battery (BT), thereby (dis-) charging the buffer E b .Due to the two clutches, both CE and MG can be disconnected and stopped, thereby eliminating their friction losses.In the following sections, the models for this topology are given, with typical model parameters denoted in Table 3.  P e0 (ω e ) W friction (cf. Figure 3) P e P e (ω e ) W minimum power (cf. Figure 3) P e P e (ω e ) W maximum power (cf. Figure 3) I e 3 kg m 2 moment of inertia MG η m 0.95 -mechanical efficiency P m0 P m0 (ω m ) W friction (cf. Figure 3) P m P m (ω m ) W minimum power (cf. Figure 3) P m P m (ω m ) W maximum power (cf. Figure 3)

Internal Combustion Engine
The internal combustion engine (CE, or 'engine') is modeled as an affine relation between the fuel P f and the power output P e , often referred to as a Willans approximation [29]: with η e the indicated efficiency (see as shown in Figure 3.The CE has P e (ω e ) ≤ P e0 (ω e ), meaning that additional engine braking can be applied on top of the nominal friction, which is a feature typically available on heavy duty commercial vehicles.For P e ≤ P e ≤ P e0 , no fuel is consumed (P f = 0).The total fuel consumption E f is the integral of P f : Ėf = P f (t). (3)

Motor Generator
The mechanical output of the motor generator (MG) is modeled as: with P m the mechanical output power, P m0 the friction, P b the electrical power from the battery and η m the constant efficiency of the MG; see Table 3.The power output is limited by and is shown in Figure 3, together with P m0 (ω m ).

Battery
The battery (BT) is modeled as an integrator, with quadratic losses, see [15]: with E b the energy in the battery, P bi the internal battery power, P b the power at the terminals and β the loss constant.Note that P b > 0 discharges the battery.The effective size of the battery is limited, such that In [8], the validation of the compression-ignition engine (CE), motor generator (MG) and the battery (BT) model is described.

Mode Selection
The topology has two clutches, CL 1 for connecting the CE to the driveline and CL 2 for connecting the MG to the driveline.When the clutch is open, the respective component is disconnected from the driveline and stopped to eliminate the friction losses in the component.The two clutches create four modes M ∈ {M c , M e , M m , M o }, representing the driveline states, as defined in Table 4.If connected, then the rotational speed of the CE(ω e ), respectively, the MG (ω m ), is equal to the rotational speed ω at the gearbox input shaft, and zero otherwise.The gearbox input power P p is the sum of the connected components.The mode is controlled with I M ∈ {i c , i e , i m , i o } for the respective four modes.Mode switch dynamics are defined by the state machine in Figure 4.When a mode switch (∆M) is performed, the driveline is open (M o ) for a duration of ∆t M , during which no traction is available (P p (M o ) = 0).During the mode switch, a series of events, involving (de-)coupling and synchronization of rotating masses, cause energy losses, represented by the lumped parameter M : In Appendix C, the dependency of M on the driveline state and the drive cycle is described.

Gear Selection
The gearbox with final drive and wheels (GB) is modeled as a geometrically stepped transmission without losses.The ratio of the gearbox r G is a function of the gear position G ∈ {1, 2, ..., 12} with the The gearbox input shaft speed ω [rpm] is related to the vehicle speed v [m/s], through the gearbox ratio, final drive and wheels, with where r d is the drive ratio from speed [m/s] to gearbox out [rpm].The gear selection is controlled with I G ∈ {i 1 , i 2 , ..., i 12 }.Gear switch dynamics are defined by the state machine in Figure 5.When a gear switch (∆G) is performed, the driveline is open (G o ) for a duration of ∆t G , during which no traction is available: During the gear switch, a series of events, involving (de-)coupling and synchronization of rotating masses, cause energy losses, represented by the lumped parameter G : In Appendix C, the dependency of G on the driveline state and the drive cycle is described.

Vehicle
The dynamics of the vehicle (VH) are modeled for the longitudinal motion: with vehicle mass m, vehicle position s, v(t) = ds dt the vehicle speed, and the total vehicle road load F 0 : with: respectively, the air drag, rolling resistance and gravitational force, and F d the driveline force: where α is the road inclination, and the other parameters as defined in Table 3.The velocity of the vehicle is limited, such that 0 with v the maximum speed limit.

Multi-Level Energy Management
Optimizing the energy consumption for a system, with multiple states and a mix of continuous and discrete decision variables, is computationally demanding for real-time implementation.Partitioning the optimization into smaller problems reduces the computational burden.Not only are the sub-problems smaller and thereby easier to solve, but also each partition can have its own optimization algorithm, making the selection of more efficient algorithms possible, suited to the specific problem of that partition.e.g., in [19] Dynamic Programming is used for the partition with discrete decision variables, while convex optimization is used for the partition with continuous states.The method of partitioning is, however, not unique nor trivial [27].
For the generic energy management problem, described in Section 3.1, a functional hierarchy is introduced in Section 3.2, inspired by the 'multi-level structure using conjugate variables' in [27].In this section, the partitioning in levels, and the solution of the multi-level EMS, for two solution schemes, are described: • the iterative scheme (Section 3.3), used in simulation to show convergence of the EMS in Section 7.1, • the model predictive scheme (Section 3.4), that is real-time implementable and simulated in Section 7.2.
The control system for each individual level is described in Sections 4-6.

Generic Energy Management Problem
The general task of the Energy Management System (EMS) is to minimize the fuel energy E f needed to move the vehicle from distance s 0 to s f .For the vehicle model in Section 2, this EMS can be posed as a nonlinear, mixed-integer, input-and state-constrained dynamical optimization problem: with continuous decision variables P d (driveline power), P b (battery power) and discrete decision variables continuous states s (distance), v (velocity), E b (battery energy) and discrete states M ∈ {M c , M e , M m , M o } (mode), G ∈ {1, 2, ..., 12} (gear): disturbances, representing the preview information, α (road slope), v (maximum velocity): switching cost M (for mode switch), G (for gear switch): equality constraints on s 0 (initial distance), s f (final distance), E b (battery energy) sustenance: inequality constraints on v, E b upper and lower limits, P e (CE power) upper and lower limits, P m (MG power) upper and lower limits: together with the model equations in Section 2.

Multi-Level Optimization
In the multi-level EMS, the global optimization problem is subdivided into three levels based on its function: velocity determination, optimizing battery energy and optimizing the power-split including switching of modes and gears.Each level (indicated with subscript i ) has the same objective, i.e., minimizing (equivalent) fuel: but with a subset of the decision variables, and on each level a different model complexity, belonging to the abstraction and dynamics on that level.In Table 5, the sub-problems of the optimization are defined, including the output y of each level.Numbering of the levels start at 1 for the lowest level, as control of the component itself (e.g., for safe operation) is indicated with Level 0, and is not contained in the EMS.
Information exchange between the three levels is shown in Figure 1 and the function and interfacing of each level is briefly explained next.
Table 5. Generic optimization problem, subdivided into three levels (i).

'Generic'
'Velocity' 'Battery Energy' 'Power-Split' Level 1 is the lowest control level in the EMS, responsible for the optimal power split (defined by P b and given P d ), mode (I M ) and gear (I G ) selection.Based on information from the higher levels, it has to act fast, in order to have a responsive vehicle.However, the required fast update rate (typically 10-100 Hz) limits the computational time available for calculations.This is solved by moving computational expensive calculations to higher layers, and having explicit solutions for the remaining optimizations, using optimal control techniques (DP and PMP, see Appendix A).Information from the higher levels is provided by vectors over time (indicated in bold) with estimated quantities: • λ b , estimated equivalent cost of battery energy, • P d , estimated power demand, • v, estimated vehicle speed, and, after optimizing the decision variables, the calculated driveline speed ω is provided to Levels 2 and 3.

Level 2, 'Battery Energy' (Section 5)
Level 2 is responsible for optimizing the battery energy over the cycle, considering the battery constraints.As the battery energy dynamics are slower than the decisions needed in Level 1, its optimizations can run at a lower rate (typically 1 Hz), facilitating the more complex calculations, caused by battery limits and longer horizons.Abstraction of the Level 1 model, by e.g., not considering gear switching and no penalties on switching, combined with PMP techniques, make the optimization computationally efficient.Based on an estimated driveline speed ω, the estimated power demand P d and the current battery energy E b , the optimal λ b is calculated, which is input for Level 1.

Level 3, 'Velocity' (Section 6)
Level 3 determines the velocity v.Of the three levels, it uses the most abstracted model of the hybrid driveline, by e.g., not considering battery dynamics, hybrid modes, power-split or switching.Using-distance based route information for preview, i.e., slope α(s) and maximum speed limits v(s), and current vehicle velocity v and position s, an estimate speed v is determined, resulting in an estimated power demand P d .

Multi-Level Iteration
The multilevel iteration scheme is used for offline simulation.In this multilevel, nested, optimization, the decisions on one level impact the objective on other levels and thereby influencing optimality [21,30].This dependency can be seen in Figure 1, where the variables and ω, λ b are fed back to the higher level controllers.Solving the optimization by starting at the highest level and sequentially transmitting information to the lower levels, we need an a priori estimate of ω and λ b .When all levels are calculated, the sequence can be repeated, where ω and λ b are updated from the last sequence, thereby providing the higher levels with the latest decision details of the lower levels.The solution scheme of the optimization results in: 1. Initialization: (a) all levels optimize over the complete drive cycle, from t(s 0 ) to final time initialize ω with an estimated average speed,

2.
Velocity prediction: (a) retrieve α(s), v(s) (road slope and speed limit information), s(t 0 ), v(t 0 ) (position and speed) and model parameters (Table 3), (b) use the vehicle model and results from optimal control, to predict the velocity profile, (c) store resulting P d , v.

Battery energy optimization:
(a) retrieve P d , ω, E b (t 0 ) and model parameters (Table 3), (b) calculate PMP necessary conditions and optimize the remaining boundary value problem(s), (c) store resulting λ b .

4.
Power-split and switch optimization: (a) retrieve λ b , P d , v and model parameters (Table 3), (b) calculate PMP necessary conditions and solve the remaining DP problem, (c) store resulting ω and output P d , P b , I M , I G to the component controllers on the vehicle.

5.
Iterate: (a) repeat from step 2 with an improved estimate of ω and λ b , unless stop conditions are met (number of iterations).Note that the initial conditions at t 0 remain identical.
It is expected that each iteration will increase fuel efficiency, until the solution is converged.The effect of the number of iterations on the fuel efficiency, is analyzed by simulation in Section 7.

Multi-Level Model Predictive Control
To apply the EMS in real time, the multilevel optimization is implemented as a model predictive (or receding horizon) controller (MPC), which performs the optimization in a regular schedule ('sampling').Each sampling period the optimization is performed with updated states, over a subset of the cycle ('horizon'), thereby creating feedback for disturbance rejection and model mismatch compensation.In Figure 6, the differences with the multi-level iteration, in horizon and iteration, are illustrated.Compared to the solution scheme in Section 3.3, the following adjustments are made: • t 0 now refers to the vehicle's current time and t f is relative to t 0 , with a fixed horizon length (which implements the receding horizon), • each sample time, only one iteration is applied, i.e., step 5 in the iteration scheme is skipped, • the horizon length decreases on each level: each level runs in its own regular schedule (sample time), where the higher levels run slower than the lower levels.
All adjustments are for improving computational efficiency; however, optimality will reduce.In Section 7, we will show that, in simulation, the optimality is marginally decreased and the multi-level approach shows good fuel economy in a high fidelity simulation environment.In Sections 4-6, the algorithms are described in detail for Levels 1, 2 and 3, respectively.The multi-level iteration has an identical horizon, for all levels, and for all iterations, comprising the complete drive cycle.The multi-level MPC has a solution scheme, with shorter horizon lengths.With each time sample, the vehicle moves along the drive cycle.At each time sample, one iteration is performed, over a horizon starting at t 0 , having a fixed horizon length to t f .When the level has a faster sample time, iterations are performed more often, thereby improving the feedback performance.

Level 1: Power-Split and Switch Optimization
The Level 1 optimization calculates the optimal power-split, mode and gear selection for a simplified plant model, as a subset from optimization (21).In [8,15], it is shown that the power-split problem with battery dynamics can be efficiently solved using PMP, resulting in two solutions steps, see Appendix A: • minimization of the Hamiltonian as a function of λ b , resulting in the optimal power-split, mode and gear, • calculation of λ b that complies to the battery constraints.
Level 1 performs the first step: minimization of the Hamiltonian.The second step, calculation of λ b , is performed at Level 2, as described in Section 5.
On Level 1, we assume a predetermined velocity profile v(t) = v, which defines P d using α(s) in Equations ( 14)- (19).As λ b is controlling E b (see Section 5), all continuous states (E b , s, v) are removed from Equation ( 23), together with the equality and inequality constraints on the respective states in Equation (26) and Equation (27).As a result, the minimization Equation ( 28) for Level 1, denoted with subscript 1 , is with with the switching cost P 1 = 0 in Sections 4.1 and 4.2 and with P 1 ≥ 0 in Section 4.3.

Explicit Minimization of the Hamiltonian per Mode
On Level 1, the Hamiltonian is minimized.For the hybrid drive train, this Hamiltonian is solved explicitly in [8,15], i.e., an analytical expression is found for the minimization.This section extends that solution for the driveline topology with an additional clutch between MG and GB.The Hamiltonian of Level 1 to minimize is using the PMP conditions (Appendix A).When we assume instantaneous switching, i.e., ∆t M = 0, ∆t G = 0, the control signals I M and I G bring the system immediately to its corresponding state, hybrid mode M, respectively, gear G, such that I M ≡ M, I G ≡ G.Then, for each M, G, the minimum of H 1 is found by solving dH 1 dP b = 0.For each M, G, the result is given in Table 6, where only parameter P m0 (ω(G, v)) depends on the rotational speed ω, as defined by G and v through Equation (11).engine-only 0 With the explicit expression for P * b (where the superscript * denotes the optimal solution), the reduced Hamiltonian H 1a (with subscript a denoting the variant) is an explicit expression, for each mode and gear combination.Finding the optimal H 1 is thereby reduced to: In Figure 7, H * 1 is illustrated.As a function of λ b , H * 1 changes modes at B 1 , B 2 and B 3 .Due to the mode changes, the control signal P * b for M * jumps as a function of λ b at B 1 , B 2 and B 3 , as shown in Figure 8.That jumping behaviour, caused by the additional clutch between MG and GB, is new to [8,15] and adds complications to the controllability of E b , as will be shown in Section 5.2.

Mode Selection, without Cost on Switching
Instead of calculating all H 1a (M, G) and selecting the minimal one in Equation ( 36), a computationally more efficient approach can be taken where the optimal mode is calculated beforehand.Figure 7 shows that the optimal mode is changed at B 1 , B 2 and B 3 where the Hamiltonians of two modes are equal.Equating the Hamiltonians, using Equations ( 1), ( 4), ( 7), (36) and Table 6, results in five expressions, representing all the switching lines (guards B) where the optimal mode M * is changing.Two guards are a function of λ b ∀P d : B 2 : and three guards are a function of (P d , λ b ): B 4 : B 5 : All optimal modes M * are explicitly defined with Equations ( 37)-( 41), as a function of (P d , λ b ) and the model parameters (Table 3).In Figure 9, the guards and optimal modes are illustrated.A further reduced Hamiltonian H 1b (λ b , P d , v) is now obtained in three steps: 1.
for each (feasible) gear, calculate select G * that minimizes H 1a (G), thereby reducing Equation (36) to: which statically defines the decision variables P b , I M , I G as a function of λ b , P d , v.
Note that power constraints of the components in Equation ( 34) are easily added by including the corresponding guards, thereby refining the area of M * , see [15].The cost function in Equation ( 29) is now reduced to which is the integral of all locally minimized Hamiltonians.

Mode Selection, with Cost on Switching
When the system operates near the guards, frequent switching (so called 'hunting') between modes and gears can occur under influence of disturbances, thus preventing acceptable real-life implementation.Each switch involves connecting and/or disconnecting of components, synchronization of their speeds and additional friction losses, which cause driveability, durability and efficiency issues.As a solution to hunting, the cost function in Equation 21 adds a penalty P = M + G ≥ 0 on mode switching ∆M and gear switching ∆G.
The parameters M , G can be tuned to balance hunting behavior, with fuel economy.To prevent tuning of the parameters, Appendix C quantifies the model-based, equivalent fuel losses during a switch, identified by: which is the integral of the explicit Hamiltonian solution in Equation (36), with switching costs M and G added to the integral, when a switch occurs.By discretizing time, with the discrete time vector .., t f 1 } and sampling time T k , the optimization is solved using Dynamic Programming.The outline of the algorithm is as follows, as illustrated in Figure 10: 1. all feasible modes M and gears G over the time horizon t f 1 are enumerated, 2.
for each combination of mode and gear and time, the Hamiltonian is explicitly solved using Equation (36), resulting in H(M, G, k), 3.
a Dynamic Program (DP) is formulated with H(M, G, k) as elements for the cost-to-go matrix (M, G, k), 4.
to implement the cost for switching, for each mode change and gear change, the cost-to-go matrix is penalized with M and G , 5. the optimal sequence of modes and gears is calculated using DP, 6.
for the MPC implementation, only the first control action is implemented.The next sample time, the algorithm is repeated with updated inputs, disturbances and states.
For the resulting M and G, P * b is selected using Equation ( 36) and with P d given, P e and P b are known with Equations ( 4), (12) and Table 4.With the selection of T k such that ∆t M ≈ ∆t G ≈ T k , we have I M (k) = M(k + 1) and I G (k) = G(k + 1), which completes the calculation of the output of Level 1.

Level 2: Battery Energy Optimization
The Level 2 optimization determines λ b , which is an input for Level 1 control.That λ b is found by solving the power-split problem, taking the battery dynamics E b in Equation (23), with corresponding limits (E b , E b , in Equation ( 27) ), into account.
On Level 2, we assume: • a predetermined velocity profile is available, which defines P d using Equations ( 14)-( 19).• a prescribed rotational speed ω, which defines G ∈ [1,12] with the predetermined velocity profile.Before iteration over the levels an estimated average ω is selected.After each iteration, ω from Level 1 is used.• mode and gear switching is instantaneous (∆t M = 0, ∆t G = 0), without costs associated (P = 0).
As a result, minimization of Equation ( 28) for Level 2 is with Following [31], this dynamic state constrained problem is solved using PMP (Appendix A).Using Equations ( 6) and (7), the system dynamics are: and with Equation (3) result in the Hamiltonian: where λ b is the costate of E b .As no switching costs are defined, H 2 is explicitly minimized with Equation (42) as a function of λ b , among others.The dynamics of λ b in Equation (A6), show that λ b is constant for the problem without constraint in Equation ( 8), or changes stepwise where E b is constrained [8,32].Solving the problem is reduced to a boundary value problem [1], with λ b the decision variable.The boundary value problem, for the topology with clutches, state constraints and a limited horizon, has the following properties, which are handled by the developed algorithms in this section, with the flowchart in Figure 11: I1 Constraint activation increases the number of boundary value problems to solve, where the number of sub-problems is not known a priori.The iterative solution method from [32] is recapitulated in Section 5.1.I2 Controllability of E b (t f 2 ) is not a continuous function of λ b , due to the clutches in the topology, resulting in non-unique solutions.The switching algorithm in Section 5.2 ensures controllability for all λ b .I3 When the horizon does not include t f , the end condition on E b (t f 2 ) is not known, thereby removing charge sustenance ∆E b = 0 as a control target.In Section 5.3, a solution is proposed, which is guaranteed to be optimal for a subset of scenarios.
In the remainder of this section, t f 2 = t f for ease of notation, but without loss of validity of the Level 2 problem where t f 2 < t f .For dealing with constraints, non-unique solutions and a limited horizon, the flowchart shows the procedure for solving the boundary value problems (BVP).Starting with one unsolved (sub)problem (SP), from t 0 to t f , the process is started.Dependent on the number of constraint violations, the procedure ends with a set of solved subproblems, each providing a segment of λ b over the complete horizon.Each subproblem contains the solution to a BVP, and dependent on the stated conditions, the additional procedures I1, I2 or I3 are needed.

Constraint Handling (I1)
The optimal solution to Equations ( 21)-( 27) results in a constant λ b when E b constraints are not present in Equation (54), which can be efficiently found by solving the two-point boundary value problem.When constraints on E b are violated, the original problem must be subdivided into one or more subproblems dependent on the amount of constraint violations, as described in [32].As we elaborate on this procedure in the next sections, we repeat an outline of the procedure here:

•
The unconstrained subproblem is solved for a constant λ b .

•
If no constraint violation is present, the subproblem is solved.The procedure is repeated for the next unsolved subproblem.

•
If constraint violations are present, the time t c of the largest constraint violation is identified.The subproblem is divided at t c into two new subproblems, with the first subproblem, having E b (t c ) = E b as an end point constraint and the second subproblem having This sequence is repeated, until all subproblems are solved, resulting in a piecewise constant λ * b for the whole cycle.

Controllability of E b : Non-Uniqueness (I2)
The optimal solution has a non-continuous E b (t f ) controllability as a function of λ b , leading to problems in realizing charge sustenance.The cause are the guards in the optimal solution of H (see Section 4.2), where a switch from one mode to another is enforced.When P d and λ b coincide with one of the guards B, the optimal solution is non-unique (singular), i.e., several control inputs lead to the same fuel optimum [31,33,34].
The first type of guards are a function of P d (and λ b ), i.e., B 3 , B 4 , B 5 (see Section 4).As P d normally varies in time, the optimal solution will seldom be at one of the guards for a prolonged duration.For that reason, Ref. [15] selects one mode a priori, when P d is on the guard, thereby being optimal, but without a guarantee to be charge sustaining.In simulations, especially with perfectly constant P d , switching between the modes could be necessary for a charge sustaining solution.Therefore, Ref. [18] opts to vary λ b around λ * b , which causes switching between the modes, but is inherently sub-optimal because of the deviation from λ * b .In [33], this is correctly solved by using a so-called 'sliding mode' control, which minimizes the number of mode changes, which does not change λ * b .In our EMS, we also unalter λ * b and remain in the last mode, when P d coincides with the guards B 3 , B 4 , B 5 .The second type of guards is not a function of P d , only of λ b , i.e., B 1 and B 2 .As λ * b is (piecewise) constant, the solution can be on B 1 and B 2 for a prolonged period of time.Figure 8 shows that the optimal solution switches from M c to M e when λ b is at guards B 1 or B 2 .Consequently P * b jumps at B 1 and B 2 , and, for an arbitrary drive cycle, the controllability of E b (t f ) (or charge sustenance) jumps as a function of λ b -see Figure 12.To have continuous controllability of E b (t f ) in B 1 and B 2 , a switching sequence must be chosen and, in the following sections, algorithms for the two switching scenarios are explained: I2a with an infinite number of switchings, I2b with a minimal number of switchings (elaborating on [33]).
The switching scenarios are illustrated for a drive cycle of 2000 s, see Figure 13, for which λ * b coincides with B 1 .Trajectory B + 1 (dash-dotted line) starts at E b (0) and is above charge sustaining at t f = 2000 s.Trajectory B − 1 (dashed line) starts at E b (0) and is below charge sustaining at t f = 2000 s.In the top plot, P m is shown for both trajectories, where the difference is in supporting the CE with MG or not (boosting).The charge sustaining solution consists of switching between B + 1 and B − 1 over time, such that E b (0) = E b (t f ).At least one switch must be made, i.e., at t s1 or t s2 dependent on the initial mode.The next subsections provide algorithms to determine the switching sequence, taking constraints on E b into account.The algorithms, as illustrated for B 1 , also hold for B 2 . (P m is zero during high power demands, resulting in too high E b (t f )).A charge sustaining trajectory (E b (t f ) = E b (0)), is possible with one switch at t s1 or t s2 , or with infinitely many switchings indicated as B ∞ 1 .

Infinite Number of Switchings (I2a)
As the dynamics of E b in Equation ( 6) are state independent, a charge sustaining solution can be realized by taking a linear combination of the two trajectories B − 1 and B + 1 , such that: where r is a ratio between 0 and 1.In the example in Figure 13, this combined solution is charge sustaining with r = 0.46 (dotted line).This combined solution implies an infinite amount of switchings between the two trajectories, and will be referred to as B ∞ 1 .

Minimal Switching (I2b)
To enforce minimal switching, considering constraints on E b , a combination of algorithms I1 and I2a is needed.First, it is checked if the subproblem must be further subdivided as in Section 5.1 to guarantee E b within its bounds: b exceeds limits, then procedure I1 for constraint handling is performed, • else a switching sequence exists in order to maintain E b within bounds.
The switching sequence with minimal switchings is determined, by choosing when to switch between the two state trajectories B − 1 and B + 1 .The algorithm is illustrated in Figures 14 and 15.

Limited Horizon (I3)
When E b (t f ) is prescribed, the boundary value problem is fully defined and the optimal solution can be found.However, in real life, E b (t f ) is only prescribed if the complete drive cycle is used (full horizon) and E b (t f ) has a value, e.g.,

•
for a repeating cycle, with E b (t 0 ) = E b (t f ), which is often used in simulation assessments, • when entering a zero-emission zone at t f , with E b (t f ) = E b for maximal electrical driving after t f , • when arriving at a charging station at t f , with E b (t f ) = E b for cheap charging of the battery after t f .
If the prediction does not include t f , e.g., when implemented as a receding horizon algorithm, The following, suboptimal, procedure aims at finding a constant λ b , without constraint violations, resulting in a E b trajectory that has maximal robustness towards constraint violations.Here, it is chosen to maximize the minimal distance of the E b trajectory to the constraints, in order to allow for disturbances in the prediction, without violating constraints.The outline of the procedure is: The iteration stops when the increment of the minimal distance is smaller than a tolerance value.When the resulting E b violates constraints, the problem is subdivided as described in Section 5.2 and each sub-problem is evaluated again.
For two special cases, the above procedure results in optimal solutions: I3a if within the horizon, a E b limit is activated, irrespective of the control action e.g., due to a large energy recuperation event, the first subproblem up to this event is fully defined.The solution to this first subproblem is then optimal, and independent from the following subproblems.I3b if (λ * b , P d ) is on guard resulting in non-unique solutions and the non-unique solutions are able to span With these algorithms, the calculation of the output of Level 2 (λ b ) is complete.

Level 3: Velocity Prediction
The Level 3 functionality determines the velocity profile of the vehicle, given the slope and velocity limitations of the drive cycle, and provides the Level 2 and Level 1 algorithms with an estimated power demand P d and velocity v.This section describes a simplified solution using results from the optimal control formulation.

Velocity Prediction Using Three Driving Modes
The Level 3 functionality predicts P d from information of the road ahead, i.e., road slope α and velocity limitations v.The velocity of the vehicle determines the resulting road load, so assumptions have to be made how the vehicle will accelerate/decelerate over the route.To determine the fuel optimal velocity profile, an optimization has to be solved.In [35], this optimization is performed for a conventional driveline with an affine fuel map, resulting in a set of three possible driving modes: 1.
maximum acceleration with F d , 2.
minimum acceleration with F d , 3.
constant speed v 0 , with The acceleration limits depend on the capabilities of the drive-line force F d , for a given road load F 0 .When we assume a continuously variable gearbox ratio, controlled to a fixed engine speed, and only one hybrid mode M e , Equations ( 2), ( 11), ( 15)- (19) and Table 4, define the acceleration limits as a function of velocity and road slope: F d (v, α) and F d (v, α).For simulation, the vehicle model is completed with Equation ( 14).We use the three modes, to predict the velocity profile, given the road information and the vehicle model.The outline of the algorithm is given: • retrieve α(s) and v(s) from the preview data source, for the drive cycle, • decide on v 0 , typically the reference speed of the speed control on the vehicle, • simulate the vehicle model forwards over the drive cycle, using mode 1 (maximum acceleration).
Limit the velocity to v 0 and v(s).Store the resulting vector of the forward simulation v f (s), With this procedure, preview information of road slope α and velocity limitations v is converted to a prediction of the power demand P d and velocity v for the lower levels.
A suggestion for future research is to calculate on this level the optimal velocity profile, considering hybrid driving modes, braking, discrete gear shifting and non-affine fuel maps.

Simulation Results of the Multi-Level EMS
Fuel optimality of the multi-level EMS is shown in this section, by three simulation scenarios: • 'Multi-level iteration, short cycle' (Section 7.1), using full horizon optimization on a short, abstract drive cycle.The optimization iterates over the levels, until the fuel consumption is converged to an optimum.With the full horizon, each level has complete information of the drive cycle.
The results show how fast the algorithm converges, and, by using an abstract cycle, the decisions on each level are interpreted.

•
'Multi-level MPC, short cycle' (Section 7.2), using receding horizon optimization on a short, abstract drive cycle.By limiting the horizon and implementing the iterations over a receding horizon, the optimality will reduce.This section shows the fuel penalty of this approach, and by using the same abstract cycle, the decisions on each level are compared to the previous scenario.

•
'Multi-level MPC, long cycle' (Section 7.3).To show the real-life fuel saving potential of the multi-level MPC, a representative, long-haul, drive cycle is used, replacing the abstracted cycles of the previous scenarios.Furthermore, a high fidelity model of the vehicle replaces the model from Section 2. The results are compared with a baseline EMS, which uses heuristics and no preview.

Multi-Level Iteration, Short Cycle
The first two scenarios are demonstrated for a cycle of 1055 m, with a vehicle weight of 25 tons.The elevation profile of this cycle is shown in Figure 16, and contains four slopes of respectively +6%, −12%, +12% and -6%.This elevation profile is also available as a physical test track as described in [28].The maximum feasible speed on this track is 25 km/h.The plant model uses the models from Section 2, including the effect of open driveline during switching (∆t M > 0, ∆t G > 0).
The multi-level optimization is iteratively solved for the complete test cycle.Multiple iterations are performed; however, no decrease in fuel consumption is observed by consecutive iterations.The energy consumption of the first two iterations are shown in Table 7, where the energy consumption difference (0.1%) is negligible small.The control decisions of these two iterations are analyzed and compared in Appendix D and show that the calculated decisions after iteration 0, hardly change with iteration 1, with the difference in charge sustenance explaining the marginal change in energy consumption.

Multi-Level MPC, Short Cycle
For real-time applications, the Dynamic Program in Level 1 over the complete route is computationally too demanding.In this section, we simulate Level 1 as receding horizon optimization.Here, we choose to optimize Level 1 at a sample time of 1 s, with a horizon of five samples.The other levels are unaltered, as they don't pose computational problems over the cycle's horizon.
The energy consumption is marginally different (0.2%) from the other scenarios.The control decisions of the multi-level MPC are analyzed and compared with the multi-level iteration in Appendix E, showing equal decisions as the multi-level iteration.One essential difference in gear selection can be explained from the shorter horizon, however, without a relevant impact on energy consumption.Comparing the three different scenarios, we conclude that the multi-level iteration scheme converges very fast: an update of the estimate of ω is not needed to improve the results of the optimizations on the layers, as iteration 0 is already at a minimal energy consumption.The multi-level MPC has a shorter horizon on Level 1, which leads to slightly different control decisions, but does not significantly alter the energy optimum.For fuel evaluation, the next section presents more realistic long-haul simulations, with the multi-level MPC implemented.

Multi-Level MPC, Long Cycle
To show the fuel saving potential of the multi-level EMS, a typical long-haul route is simulated, in which a non-hybrid vehicle is compared to a hybrid vehicle with baseline EMS ('baseline') and the multi-level MPC ('preview').The baseline EMS is a proprietary non-previewing, heuristic algorithm.The proprietary simulation environment is designed for fuel evaluation purposes and contains high-fidelity models of the vehicle.The multi-level MPC is implemented, with the following settings: For the 25-ton vehicle, the fuel consumption results are shown in Table 8.The difference between 'baseline' and 'preview' is best illustrated with the fuel saving and battery energy E b , plotted in Figure 17.There we observe the largest advantage of the previewing strategy, when battery limits are frequently touched, i.e., in the hilly part of the cycle: from 150-350 km.To recuperate a maximum amount of brake energy, 'baseline' discharges the battery as soon as possible, which is not the most efficient use of the recuperated energy.The previewing strategy discharges the battery sufficiently, to recuperate the maximum amount of brake energy, and uses the stored energy more efficiently, e.g., to drive MG-only during low road loads.Also on segments without battery limitations, the previewing strategy outperforms the baseline strategy.On the flat road, e.g., from 0-120 km, λ b coincides with guard B 2 (see Figure 9 and algorithm I3b in Section 5), resulting in fuel optimal switching control: charging the battery is performed in a relatively short period of time, shown by the steep increase of battery energy, followed by a prolonged period of a disconnected MG, thus eliminating MG friction P m0 .The baseline strategy, however, uses lower charging powers, thereby having little periods of MG disconnection, and thus higher MG friction.Furthermore, from a driveability perspective, the larger duty cycles of the previewing strategy, in connecting/disconnecting the MG, are preferred.

Conclusions
In this work, a computationally efficient EMS using preview is developed, obtaining excellent fuel economy results (−7.0%) on a typical long-haul route.To manage computational complexity, a multi-level EMS is designed, with on each level a specific algorithm for the subproblem to solve.Each subproblem uses different model information of the vehicle, such that optimal control techniques can be used, to efficiently solve the problem.
On level 1, the power-split, gear, mode decisions are optimized, including costs on switching.For the topology with two clutches, a novel explicit PMP solution is found, which reduces power-split and mode selection to a static calculation, dependent on the current vehicle state and battery costate (or Lagrange multiplier).For taking switching costs into account, a Dynamic Program is designed, using the explicit PMP solution for reducing the dimension of the problem, and thereby limiting the complexity of the calculations.
On level 2, the use of battery energy is optimized, taking the state constraints of the battery into account.A simpler, more abstracted vehicle model is used, e.g., by ignoring switching costs, to reduce computational complexity.Using PMP, the battery costate (or Lagrange multiplier) is calculated, which provides essential cost information to Level 1.The topology with two clutches, causes discrete jumps in the optimal solution, leading to non-unique (singular) solutions.New optimal solutions are found, minimizing the amount of switching events.
On level 3, the velocity is predicted along the route, using road slope and velocity limitations.With the calculation, essential route information is provided to the lower layers.Future work should include velocity optimization on this level, to further reduce the fuel consumption.
In offline simulation, on a short abstract drive cycle, the multi-level approach is shown to converge within one iteration.When the multi-level approach is used as Model Predictive Controller (MPC), comparable fuel energy efficiency is realized, even though the horizon is reduced.
Representative fuel saving results are obtained by running the multi-level MPC in a high fidelity simulation environment, over a typical heavy duty drive cycle: 7.0% fuel reduction compared to a conventional vehicle, where the baseline EMS for the hybrid electric vehicle saves 5.8%.

•
each switch involves an inverse switch on a later moment.
The last assumption makes it possible to allocate an average cost to two linked events, e.g., starting the Combustion Engine (CE) has large costs associated, while stopping the CE is free.However, the one can not occur without the other, so the total cost of the two events can be averaged over the two single events.This results in balanced switch decisions in the EMS, when the horizon is small and only one of the linked events can be overseen within the horizon.
With the inertia from Table 3, the equivalent cost of the transitions are quantified with the change of rotational kinetic energy: and translated to the equivalent fuel with: which results in with the inertia and speed just before (I − , ω − ) and after (I + , ω + ) the mode switch.As λ b defines the fuel equivalence between electrical and thermal energy, it is not relevant with which energy source the synchronization is performed.For example, an ICE start first involves an electrical acceleration, followed by mechanical acceleration as soon as the rotational speed is above idle speed, but, with λ b , they are equivalent.The same method is used to calculate c Gs .
For the following rotational speeds, • disconnected, stopped: 0 rpm, • connected: 930 rpm, • downshift: 1200 rpm, the equivalent costs (here: λ b = −2.2) are exemplified in Table A1.The downshift event losses ( * ) depend on the active mode, respectively M m , M e , and M c , because losses increase with the amount of rotational mass to be accelerated, which is the largest with both CE and MG connected.The transitions involving deceleration of the rotational masses are set to 0. However, some energy recuperation could be implemented with the MG, resulting in lower average equivalent costs.It should be noted that the transitions are idealized.The real energy consumption of the transitions will be worse than the posed numbers due to additional losses but can be improved by recuperation on the decelerating masses.Nevertheless they provide a first estimate on the selection of , which has clear mode dependent values.For a more thorough analysis of the gear switch event, including clutch slip, see [41] (Chapter 6).

Appendix C.2. Traction Interruption Cost
During the transition of a switch, the traction is momentarily interrupted.This interruption causes a deviation from the intended vehicle speed, which has to be compensated during the periods with traction.Dependent on the road load, the cycle and the component capabilities, this compensation is expensive or cheap.The idea here is to quantify the cost of this compensation to the cycle through the information contained in the costates belonging to the cycle: λ t , λ k in Equation (A14).To link this information to the traction interruption cost, the following idealized transition is assumed: • the transition has a fixed traction interruption time (∆t o ), • during the transition, the road load is constant, • the driveline is not saturated, i.e., traction is available to recover the vehicle speed, • λ t represents the equivalent cost of time loss for a certain cycle, and is constant.
In Figure A1, the idealized transition is plotted.While cruising at v 0 , where F d = −F 0 , a switch event is planned at t 0 .During the following open driveline period t 0 -t 1 , constant deceleration causes the vehicle speed to drop.After the transition, we assume immediate recovery to the desired speed, with constant acceleration at F d , depending on the mode M and associated component capabilities.During this event, the vehicle loses time, with respect to cruising without transition.If we compare the total transition+recovery time to reach the position s 2 at t 2 , with the time needed for constant cruising to reach position s 2 , we define the time loss ∆t t of the transition event: with This time loss is converted to equivalent fuel, with λ t , thus providing a road load dependent penalty on switching, When the drive-line saturates (|F 0 | ≥ F d ), the vehicle is not able to recover the cruise speed after a short open driveline event.In that situation, the speed deviation v 1 − v 0 is used with λ k to calculate the speed recovery cost.The same method is used to calculate c Gt .
As an example the equivalent cost of a transition due to time loss is given in Table A2, for the model parameters from Table 3, the vehicle driving 25 km/h, and the cycle having a λ t of −25 kJ/s.The event is tabulated for three segments of, respectively, 0%, 6% and 12% slope.The traction interruption period t 1 − t 0 is set to 1 s.
Table A2 shows that more than a magnitude order of difference occurs for a switch cost.This means that the influence of time loss on the equivalent fuel consumption is relevant when the open driveline deceleration (acceleration) is high, where on other route segments the influence is negligible.This shows that transitions can be best planned during the segments with smallest deceleration (acceleration).
The switch cost is now determined by two elements: synchronization cost as a function of the intended transition, and time loss cost as a function of position along a route and the available traction for the transition, in Equations (A15) and (A16).In Figure A2, the result is shown for respectively MG power (P m ), mode (M), gear (G), battery energy (E b ) and vehicle speed (v).The top plot also includes P d (dotted black) as a reference.On the flat segments, P d = 11 kW and on the slopes P d = [113, −192, 215, −90] kW for, respectively, the +6%, −12%, +12% and −6% segment.The dashed lines indicate the decisions of Level 2, the solid lines of Level 1.  Level 2 (dashed line) calculates a charge sustaining λ b = −1.97,shown in subplot 4, where E b starts and ends at 0. It uses a simplified model which keeps the engine speed at 1100 rpm, showing in subplot 3 as a selected gear between 7 and 8.At this engine speed, it is optimal to drive the flat segments in MG mode (M m ), with the CE switched off.The uphill segments have higher power demands P d , and the CE is used.On the downhill segments, brake energy is recuperated with the MG, but the CE must assist in braking, as |P d | is larger than the MG can recuperate.As all switches are assumed to be infinitely fast, the bottom plot shows no speed deviation for Level 2.
Level 1 (solid line) uses λ b = −1.97 and the switching costs to determine the optimal mode and gear switching sequence.Most importantly, the gear on the downhill segments is shifted down to 5, in order to maximize the recuperation power of the MG (see Figure 3), which is notable by an increased rise in energy of the battery for these segments.Due to the downshift, the −6% can be driven with the MG alone.The −12% still needs both MG and CE for additional braking power.The uphill segments show a different mode selection than for Level 2. From the control map, we know that λ b = −1.97 is close to the guard G 1 , where the Hamiltonian of M e is only slightly lower than of M c , causing Level 2 to decide for M e .Level 1, however, takes the switching costs into account, and where both levels decide to start and stop the CE, Level 1 decides that stopping and starting the MG is more expensive than the minor fuel advantage of the lower Hamiltonian, resulting in mode M c with 22 kW boosting of the MG.Gear selection on the flat segments is towards the highest possible gear where the MG is able to deliver the demanded power, thereby minimizing the friction of the MG.
Due to the switching of mode and gear, several open driveline events occur.The cost of timeloss, as calculated in Section 4.3, determines that, if a switch is needed, it should be done where the timeloss is minimal.For this cycle, this is clearly on the flat segments, not on any of the slopes.This results in only minor vehicle speed deviations in the vehicle speed.
Due to the different control actions between Level 1 and Level 2, charge sustenance is not enforced anymore.This is to be corrected in the next iteration.

Appendix D.2. Iteration 1
In the first iteration, each control level is optimized consecutively, using the results from the zeroth iteration.The differences with the zeroth iteration are: • Level 3: the constant engine speed from iteration 0 is replaced by the speed resulting from the gear switching sequence.As speed is not explicitly optimized and can always be maintained in this example, P d is not altered.• Level 2: using the gear switching sequence, a new charge sustaining λ b is calculated.

•
Level 1: with the adjusted λ b , the power-split and switching sequence are optimized In Figure A3, the results are shown for Level 2 and Level 1 for the first iteration.It is evident that no significant differences are present between the two levels.An additional gear shift occurs on the top of the first hill, but, as the friction of the MG at different speeds is very small, the fuel advantage is negligible.Comparing to iteration 0, the difference with Level 2 of iteration 0 is also very small.The main effect of the iteration is that the cycle is charge sustaining again, with λ b = −1.95.To realize that, the boosting power of the MG is slightly increased from 22 kW to 29 kW on the uphill slopes.

Appendix E. Simulation Results of Multi-Level MPC
In Figure A4, the results are shown for Level 2 and Level 1 for the first iteration.As Level 2 is not changed, the Level 2 traces are the same as in Figure A2.Albeit Level 1 now has a horizon of 5 s, instead of 10 km, the mode decisions are unaltered.However, the gear decisions differ in some segments:

•
the trace starts at G 5 and is maintained for a long time.That gear is selected at the end of the −6% slope to maximize brake energy recuperation.After the slope, P d is low, and a little friction reduction of the MG could be accomplished with an upshift.However, the fuel improvement of switching in MG is an order of magnitude lower than in other modes.Now, for a gear switch to occur, the accumulated fuel improvement over the horizon (5 s) must be higher than the cost of the switch event itself, which it is not, for this short horizon, and the switch will not occur.
• on the +6%, a higher gear G 9 is needed, in order to minimize CE friction losses.Before the +6% slope, this decision enters the horizon, and, as the switch will be made anyhow, the algorithm decides to do the switch immediately, as driving MG in G 9 is slightly more profitable.

•
on the top of the first hill, a short shift to G 10 is performed, simultaneously with a mode switch, which makes the gear switch event cheap.The decision has, however, a marginal impact on the fuel consumption.

Figure 1 .
Figure 1.Multi-level optimization of the velocity (v), battery energy (E b ) and power split (P e , P b ) including mode and gear selection (I M , I G ).In Section 3 the scheme is described, with the prediction of: road slope (α), velocity limit (v), driveline rotational speed (ω), driveline power (P d ), velocity (v) and battery costate (λ b ), dependent on: vehicle position (s), velocity (v), battery energy (E b ), driveline power (P d ), mode (M) and gear (G).

Figure 2 .
Figure 2. Topology of the parallel Hybrid Electric Vehicle with its relevant power flows, defined positive towards the road and relevant states (between parenthesis).

Figure 3 .
Figure 3. Power capabilities of internal combustion engine (CE) and motor generator (MG) as a function of the rotational speed.Dotted lines indicate the nominal friction P e0 and P m0 .

Figure 4 .
Figure 4. State machine for switching between modes M. All mode switches go through M o ('open driveline'), which is maintained for ∆t M .The negation of the mode command is denoted with ∼ i.

Figure 5 .
Figure 5. State machine for switching between gears G.All gear switches go through G 0 ('open driveline'), which is maintained for ∆t G .The negation of the gear command is denoted with ∼ i.

Figure 6 .
Figure 6.Schematic of the horizon in the 'multi-level iteration' and 'multi-level Model Predictive Controller (MPC)' solution schemes.The multi-level iteration has an identical horizon, for all levels, and for all iterations, comprising the complete drive cycle.The multi-level MPC has a solution scheme, with shorter horizon lengths.With each time sample, the vehicle moves along the drive cycle.At each time sample, one iteration is performed, over a horizon starting at t 0 , having a fixed horizon length to t f .When the level has a faster sample time, iterations are performed more often, thereby improving the feedback performance.

1 Figure 7 .
Figure 7.The Hamiltonian H * 1 of the three modes M m , M c , M e as a function of λ b for a power demand and rotational speed.The optimal mode M * changes at the boundary points B 1 , B 2 , B 3 .

Figure 8 .
Figure 8. P * b as a function of λ b for the modes M m , M c , M e for a power demand and gear.At B 1 and B 2 , the optimal mode M * switches between M c and M e , with a corresponding jump in P * b .At B 3 a jump occurs due to switching between M c and M m .

Figure 9 .
Figure 9. Explicit mode change map as a function of costate λ b and power demand P d , indicating the guards B between the optimal modes M * m (motor-only), M * e (engine-only) and M * c (combined).Figures 7 and 8 are cross sections at P # d .

•
synchronization losses, caused by acceleration and deceleration of rotating masses (Appendix C.1), • traction interruption, caused by disconnection of the driveline and leading to vehicle speed deviations (Appendix C.2).As the speed deviation during disconnection of the driveline, depends on the road load (e.g., up-hill resulting in a large speed decrease, on flat road resulting in a small speed decrease), the corresponding equivalent fuel costs are a function of the speed change and time loss.The calculated penalties are then a function of mode change, gear change and vehicle speed.Adding the cost of switching, changes Equation (29) to the new Level 1 cost function min

Figure 10 .
Figure 10.Illustrating the Receding Horizon Integer Program to solve at each time step.Over a horizon of four seconds with a sampling time of one second, all feasible modes and gears are enumerated (blue), starting at the current mode and gear (red).For all feasible modes, H is calculated, increased with an equivalent fuel penalty P on each mode or gear change, and the optimal trajectory over the horizon (blue dotted line) is determined, using Dynamic Programming (DP).

Figure 11 .
Figure11.For dealing with constraints, non-unique solutions and a limited horizon, the flowchart shows the procedure for solving the boundary value problems (BVP).Starting with one unsolved (sub)problem (SP), from t 0 to t f , the process is started.Dependent on the number of constraint violations, the procedure ends with a set of solved subproblems, each providing a segment of λ b over the complete horizon.Each subproblem contains the solution to a BVP, and dependent on the stated conditions, the additional procedures I1, I2 or I3 are needed.

Figure 12 .
Figure 12.Changing the constant λ b for a drive cycle, causes jumps in the charge sustenance, due to non-unique (singular) solutions at B 1 and B 2 .

Figure 13 .
Figure 13.Illustration of algorithm I2a.Power demand P d of a drive segment with λ * b = B 1 and two optimal modes: B − 1 (P m is positive during high power demands, resulting in too low E b (t f )) and B + 1

•
choose λ b based on bi-section and simulate E b , • find the time index of E b and E b , • calculate the distance to the boundary for the two indices and select the minimal distance, • select the next λ b with bi-section, to increasing the minimal distance, and iterate.
• simulate the vehicle model backwards over the drive cycle, using mode 2 (minimum acceleration).Limit the velocity to v 0 and v(s).Store the resulting vector of the backward simulation v b (s), • calculate the minimum of the two vectors v c (s) = min(v f (s), v b (s)), • calculate P d using the vehicle model and v c (s).

Figure 16 .
Figure 16.Elevation profile of the short test cycle.

Figure 17 .
Figure 17.High fidelity simulation results on a representative long-haul cycle.The previewing Energy Management System (EMS) saves an additional 1.2% fuel to the baseline EMS.

Figure A2 .
Figure A2.Decisions and results from Level 2 and 1, after the zeroth iteration.Due to the different model abstraction, the solution on Level 1 is not exactly charge sustaining anymore.

Figure A3 .
Figure A3.Decisions and results from Level 2 and 1, after the first iteration.After updating the information on λ b and ω, the solution on both levels converges, and are charge sustaining.

Figure A4 .
Figure A4.Results for MPC, with a 5-sample horizon on Level 1. Differences with Figure A2 are mainly in the gear decisions, with a marginal impact on fuel.

Table 1 .
Optimizing Energy Management System (EMS) algorithms, for power-split, stop-start and gear.'Cost'indicatesif switching events are penalized in the algorithm.See Table2for the used abbreviations.

Table 3
) and P e0 (ω e ) the speed dependent friction losses.It should be noted that P f (t) ≥ 0 and P e0 (ω e (t)) ≤ 0. The power output is limited by P e (t) ∈ [P e (ω e ), P e (ω e )],

Table 6 .
Explicit solution of the power-split as a function of the mode.

•
The initial mode M is maintained, here resulting in B − 1 , • Check for constraint violations.E b violation is detected at t = 388s, so switching to B + 1 is mandatory.Switching at t = 388s results in E b undershoot of E ∧ b .With the known relative E b trajectories of B − 1 and B + 1 , the advancement of the switching time is determined, here at 83 s, Figure 15, • The previous step is repeated, until t f is reached, • E b (t f ) is now within limits, but not at E b (0).A switch to the other mode is back propagated from E b (t f ) = E b (0), resulting in a last switch at 1191 s.With this procedure, the discontinuity in the controllability of E b (t f ) is solved with a switching sequence between modes, where the number of switchings is minimal.The value of λ * b has not been altered.Note that knowing the switching sequence is not necessary on Level 2, as the interface from Level 2 to Level 1 does not use E b , but λ b , which is not influenced by the switching sequence.Knowing that a switching sequence exists is thereby sufficient.The Level 1 optimization then decides on switching, taking switching costs into account (Section 4.3).Illustration of algorithm I2b, first part.Starting with B − 1 the lower boundary is violated at t s and a switch must be made to B + 1 .However, t s = 388 s as switching point is too late, as B + + Figure 14.b Figure 15.Illustration of algorithm I2b, second part.A new t s1 = 83 s is found, and B + 1 is continued, until at t s2 = 1191 s the switch to B − 1 has to be made to ensure E b

Table 7 .
Energy consumption over the short test cycle.MPC: model predictive controller.

Table 8 .
Fuel savings on a typical long-haul route.EMS: energy management system.

Table A1 .
Equivalent cost of transition: synchronisation.

Table A2 .
Open driveline behavior for timeloss (∆t t ) calculation.From t 0 to t 1 , no traction F d is available, causing the vehicle to decelerate.From t 1 to t 2 , speed recovery takes place at maximum traction F d .Arrival at position s 2 (t 2 ) is later than without switch event (dotted line), which is the defined timeloss.Equivalent cost of time loss due to a switch event of 1 s.