Model-Predictive Energy Management for Optimizing the System Efficiency of a Fuel Cell Truck Considering Fuel Cell Degradation

1. Introduction and Motivation

In order to comply with the maximum global temperature increase of 2 K by 2050 set by the Paris Climate Agreement in 2015, the European Union (EU) has formulated greenhouse gas emission limits for the transportation sector [1]. Based on 2016 data, heavy-duty commercial vehicles are responsible for 27 % of $C{O}_2$ emissions from road transport and around 5 % of greenhouse gas emissions in the EU, according to a press release from the European Parliament [2]. This is due to the increase in freight transport since 1990, which has increased CO₂-emissions from heavy-duty vehicles by 25 %. Electric powertrain architectures offer a future-proof solution for meeting greenhouse gas emissions regulations in all possible transport sectors [3]. The transition from fossil fuels to environmentally friendly energy sources introduces new challenges that must be addressed to meet consumer and regulatory demands.

In addition to electrified powertrains, so-called electric vehicles (EVs), fuel cell electric vehicles (FCEVs) in particular offer a promising solution in terms of zero local CO₂-emissions, higher energy density and faster recharging compared to battery electric vehicles [3,4]. The reduction of greenhouse gas emissions can be achieved through the integration of a fuel cell and green hydrogen [5] within a powertrain application. In an all-electric FCEV, the fuel cell is the primary energy source [6]. With more than one energy source, a FCEV uses a combination of a fuel cell system and an energy storage system (ESS), preferably a HV traction battery or supercapacitors (SC), to compensate for the slow dynamics of the fuel cell during load changes and to enable the regenerative recovery of braking energy [6,7]. It is evident that fuel cell systems represent a viable option for achieving CO₂-neutral mobility, particularly in the context of heavy-duty vehicles [8]. In order to operate a hybrid system with multiple energy sources in an efficient way, an active energy management strategy (EMS) is required. In the case of fuel cell operation, the power distribution must be optimized to operate at the most efficient point, taken into account system constraints and durability [3,4].

1.1. Energy Management Strategies

EMS are operating strategies that attempt to efficiently meet various powertrain requirements, taking into account the optimal use of multiple energy sources. EMS are used to improve the performance of hybrid vehicles in terms of both energy consumption and component lifetime [7].

EMSs can be classified into three groups: rule-based, optimization-based and learning-based EMS. The proposed methods have their own advantages and limitations, which are discussed individually in the following section. Figure 1 provides an overview of the most frequently used optimization algorithms, classified according to their real-time capability and their potential for optimization over the prediction horizon.

To reduce the load change rate of the fuel cell system [9] used a series fuzzy control strategy. As a result, the load changing rate is reduced and limited to the range of change capacity. Fuzzy-logic (FL) control is an extension of conventional rule-based (RB) control based on if-else conditions that attempt to recreate the human decision-making process [10]. Another type of RB-EMS is frequency decoupling, which splits the power demand in the frequency domain so that the high-frequency demand is met by the battery and the low-frequency demand is met by the fuel cell [9,11]. This method has been extended by [12] with an optimization-based strategy to consider not only the dynamic requirements but also the energy efficiency in terms of fuel consumption and state of charge (SoC) management.

Rule-based strategies can be implemented with little effort, since they require neither a plant model nor training data, and their rules can operate on both directly measured and model-based signals. They have a low computation time due to the stationary control decisions based on the input values. However, the heuristic strategies require a high calibration effort in order to guarantee a performance close to optimality for each driving cycle. In addition, the rules are not easily scalable to different powertrain architectures or different component sizes [13].

Learning-based (LB) EMS use data mining methods that utilize historical and real-time information from large amounts of data to derive the optimal control law.

An intelligent energy management concept that employs training data to determine the best power control for a hybrid vehicle system with three or more energy sources and significant load fluctuations is presented in [14]. The machine learning algorithm covers a wide range of power requirements with restrictions in terms of system voltage and keeping the states of charge of the energy storage units within predefined thresholds. By optimally combining the energy sources, the power losses of the overall system were reduced by 28 % up to 76 % compared to a conventional rule based power control [14].

Reinforcement learning (RL) is presented in [15], taking into account the lifetime of the fuel cell system. The RL algorithm is pre-initialized by energy allocation rules to speed up the optimization process. Simulation results show that the average H₂ consumption difference between the proposed EMS and a EMS based on DP is 5.59 % in the training and validation driving cycles. In addition, the proposed EMS reduces average power difference of the FCS, which should also reduce the performance degradation. Compared to a conventional RL-based EMS, which does not take into account the lifetime of the fuel cell system, the power fluctuations could be reduced by at least 13 %. This is a significant advantage for improving the lifetime of the fuel cell system.

The learning-based strategies using the example of RL require a large amount of historical data for the development of control strategies and are correspondingly computationally intensive [7,16]. They provide a solution close to optimality and are particularly suitable for multicriteria optimization problems [16]. However, LB strategies based on deep neural networks, such as RL, operate as so-called black box models: the internal decision-making process of the trained network is non-transparent and cannot be formally verified or traced back to physical quantities by engineers [17]. This inherent opacity is fundamentally incompatible with automotive functional safety standards such as ISO 26262, which require formal verification and requirement traceability for safety-critical software in series production vehicle control units [18]. Consequently, LB-based EMS are currently not applicable for series production control units.

Optimization-based (OB) algorithms regulate the control variables based on numerical computation results by minimizing a predefined cost function with feasible constraints [19]. According to [7], Dynamic Programming (DP) and Pontryagin’s Minimum Principle (PMP) are the most frequently used global optimization algorithms. In [20], the DP for a FCEV is used to optimize the fuel cell power to reduce hydrogen consumption and extend the lifetime of the fuel cell. An optimal power distribution strategy for different degradation stages of fuel cells was derived in [21] using DP. Compared to a rule-based strategy, the operating costs of the fuel cell were reduced by 14.17 % and the average lifetime was increased by 8.48 %.

In [22] a health-conscious two stage predictive energy management strategy for fuel cell electric heavy-duty long-haul trucks is proposed to overcome the challenges of durability and efficiency. The multi-objective optimization problem integrates component degradation as an inherit part of the optimization process to reduce the degradation of the battery and the fuel cell taking into account the performance deterioration due to the progressive degradation of the powertrain components. To do so in the first stage the fuel cell power and SoC profiles are optimized for the entire route once at the beginning of the driving cycle. The optimized references are followed by an on-board controller, representing the second stage of the predictive EMS which is based on DP as an online EMS developed in [23].

For a FC vehicle [24] combined an optimal control using the PMP with a Markov chain model. The Markov chain prediction model was used to predict the future energy demand during the driving cycle. Another study used PMP as a global optimization method for the power distribution between the fuel cell and energy storage system in order to minimize hydrogen consumption [25]. Using a developed degradation model, a factor for limiting the fuel cell’s power fluctuations with a weight coefficient was added to the PMP in [26] to increase the lifetime of the fuel cell. In comparison to a RB strategy, simulation results show a reduction of the fuel cell degradation by 56.4 % in the China_city driving cycle respectively 57.8 % in the Highway fuel economy test (HWFET) driving cycle. A similar strategy was used in [27] to limit fuel cell power fluctuations, extended by a SoC range limitation between 40 % and 80 %. As a result, dynamics of the fuel cell output could be reduced by 38.3 % compared to a rule-based strategy.

OB-EMSs such as DP and PMP are able to find a global optimum for the given problem. The DP achieves the highest degree of optimality, but requires the highest computational effort and the longest computing time. The PMP approach can find a comparable optimum with a significantly lower computational cost. However, both strategies are mostly implemented offline due to higher computation requirements. Therefore, local optimization algorithms are primarily the focus for online implementation as shown in [16,28].

Local optimization algorithms such as the Equivalent Consumption Minimization Strategy (ECMS), as well as Model Predictive Control (MPC) and stochastic dynamic programming (SDP), can find a local optimality close to the global optimal solution and can be implemented as real-time applications on control units.

The Equivalent Consumption Minimization Strategy (ECMS) is an online adaptation of the PMP algorithm where the global optimization problem is converted into a local problem.

Fuel consumption is specified as the optimization target in [29]. The total costs are made up of the equivalent cost of hydrogen consumption and the equivalent cost of the battery’s electricity consumption. Since the driving profile is not known in advance, an equivalence factor correction function is used to adjust the SoC of the battery. Here, the actual equivalent factor is calculated by multiplying the SoC correction function by the nominal equivalence factor. In [13], a two-stage ECMS is used to minimize the hydrogen consumption while protecting the PEMFC at the same time. For this purpose in the first stage, a predictive controller is developed using the telemetry equivalent consumption minimization strategy (T-ECMS) approach, which predicts the global optimality trend of the battery state of charge and the local control reference using telemetric driving information data. In the second step, a tracking controller is used that tracks the local control reference in the current control step to follow the optimal SoC reference, taking into account the PEMFC state constraints and other physical constraints. With the advancements in the availability of real-time information [19] presents a multi-layer predictive energy management with a ECMS in the lower layer combined with information such as live traffic data to reduce the fuel consumption. The proposed algorithm leads to a reduced fuel consumption in comparison to a RB strategy for a selected real-world use-case by 7.5 %.

ECMS can be easily implemented in real-time applications but have a lower adaptability for a fixed equivalence factor (EF). Their optimality is strongly dependent on the choice of EF and does not guarantee global optimality [16]. Adaptive ECMS (A-ECMS) can be used for online EF adaptation, providing local optimal solutions in real-time [30]. The real-time adaptation increases the complexity and thus the computing effort of the ECMS.

Model Predictive Control (MPC), also called Receding Horizon Control (RHC), is a global optimization algorithm that is applied as an online method for a finite time interval into the future, called the prediction horizon.

An EMS based on a linear MPC has been developed in [31] to improve battery performance and prevent fuel cell and battery degradation. Taking into account the different operating modes of the battery and the time constraints associated with the activation and deactivation of the fuel cell, an effort is made to maintain the discharge level of the battery at the optimal performance threshold and to minimize the hydrogen consumption. In order to reduce hydrogen consumption and improve the performance of the hybrid vehicle, [32] developed an MPC as a nonlinear constraint optimization problem to predict the optimal power of the fuel cell system in conjunction with prediction models for the future demand power for an FC vehicle. By using neural networks, improved results could be achieved for the optimization variables, in comparison to a Markov chain prediction model. However, neither strategy takes into account the aging of the fuel cell, except indirectly through operational constraints.

In the broader literature, fuel cell degradation is treated at different levels of detail depending on the application and the available information: (i) via operational constraints (e.g., limiting current/power gradients, avoiding start–stop events, and preventing operation in high-voltage or fuel-starvation regimes), (ii) via surrogate degradation metrics that can be embedded into an optimization objective (e.g., dynamic stress factors or voltage-decay-rate proxies), and (iii) via explicit state-of-health (SoH) states and deterioration-aware models that adapt the power split over the lifetime. In heavy-duty applications, experimental studies highlight the relevance of operating conditions for system efficiency and durability [33]. At the system-modeling level, map-based or efficiency-based representations are often used to enable integration into larger simulation frameworks, including other application domains such as aeronautics [34,35].

The different strategies are compared in terms of their optimality, computation time, degree of implementation and real-time capability (see Figure 2).

Due to the ability to achieve a near-optimal solution to the DP with a high prediction horizon, as well as the consideration of system and control constraints, the chosen EMS for this work is the Model Predictive Control approach. To ensure real-time capability, the computational requirements of the MPC are solved using the advanced Acados solver [37,38].

1.2. Problem Definition

This work deals with the development of a model predictive control to optimize the energy efficiency of fuel cells, taking into account the lifetime of fuel cells. The described control strategy in [39,40] is used as a baseline strategy for the hybrid management as well as for the control of the fuel cell system. The contribution of this work with respect to existing literature is the formulation and real-time-capable implementation of a degradation-aware nonlinear MPC for a fuel cell heavy-duty vehicle, which (i) embeds a dynamic-operation degradation proxy directly into the OCP as an objective, (ii) jointly optimizes hydrogen consumption, system-efficiency-related terms and fuel cell aging within a consistent constraint framework, (iii) is validated in a MiL environment on both a representative heavy-duty cycle (VECTO Long Haul) and a real driving cycle, with a comparison to both a global optimal DP benchmark and a calibrated RB baseline., and (iv) can be implemented on automotive control units.

The paper is structured as follows: Section 2 describes the used simulation model for the FC HDV, Section 3 delves into the development of the MPC. Section 4 provides a detailed analysis of the potential in hydrogen savings and the extension in lifetime of the fuel cell.

2. Simulation Model

In order to generate an accurate virtual representation of the real system, a detailed plant model in a MiL simulation is used. The plant model can be divided in two main parts, one being the stack model and the other the system model, composed by the auxiliary components of the fuel cell system and the hybrid system.

Depending on the system structure, the components can differ greatly from one another and therefore data maps are used to represent the different components. On the other hand, for the stack model [41], the approach is to use physical equations and a segmentation in the three spatial direction in order to achieve a high level of accuracy and account for local phenomena.

A visualization of the segmentation in the in-plane and through-plane direction of a stack cell is shown in Figure 3. On each side of a section in the through-plane direction, a bipolar plate is present through which the reaction gases and the coolant are fed. Further inwards there is additionally a section that represents the gas diffusion layer (GDL), which can be supplemented by a micro-porous layer (MPL) in order to model different types of fuel cells. In the middle the cathode catalyst layer (CCL), the anode catalyst layer (ACL) and the membrane compose the membrane electrode assembly (MEA). For improved accuracy, the membrane and each GDL/MPL can be segmented in further parts, so that local phenomena can be modeled more precisely. The different effects taking place in the stack model can be separated in three distinct sections: mass transport, thermal dynamics and electrochemistry.

The mass transportation is separated in the model in two distinctive parts. On the one hand, the movement of gases such as oxygen, hydrogen or nitrogen is considered. The transport mechanisms of water is modeled separately, as water molecules in the gas phase can change to the liquid condition depending on the conditions present in the system and the two phases interact with each other. The motion of gas molecules can primarily be described by Knudsen diffusion whereas liquid water is being transported through capillary pressure [42]. Moreover the transport of water through the membrane is being considered. The distinction between sorption and desorption and the consideration of transient phenomena are particularly important in order to be able to better model the membrane behavior [43,44].

The thermal model considers the generated heat of the chemical reaction taking place at the cathode and the heat being generated through the transport of hydrogen ions being transported through the membrane and the electrons moving through the different cell components. The heat transfer in between the two bipolar plates is taken over by heat conduction, whereas convection is the primal force of the transfer from the gas channel material to the coolant. All effects of heat generation and transport are regrouped with the heat capacities of the different materials by a system of differential equations, which are solved to determine the temperature at various positions and at different points in time [42].

The calculated variables are used as an input for the electrochemistry model. Here different effects are grouped together to determine the correct cell voltage. The theoretical cell voltage $U_{\mathrm{th}}$ can be calculated with the following equation:

$$\begin{array}{cc}\hfill U_{\mathrm{th}}& =-{\displaystyle \frac{\Delta G^0}{z·F}}-{\displaystyle \frac{R_{\mathrm{m}}·T}{z·F}}ln\prod {\left({\displaystyle \frac{p_{\mathrm{i}}}{p_0}}\right)}^{{\upsilon }_{\mathrm{i}}}\hfill \end{array}$$

(1)

The cell voltage is dependent on the Gibbs free energy of reaction $\Delta G^0$, the stoichiometric charge number z, the Faraday constant F, the universal gas constant $R_{\mathrm{m}}$, the temperature T, the pressure of $H_2$ and $O_2$ in the gas channel $p_{\mathrm{i}}$, the reference pressure $p_0$ and the stoichiometric coefficient ${\upsilon }_i$ of $H_2$ and $O_2$.

Additionally a constant voltage loss can be subtracted to obtain the OCV voltage. The other representative voltage losses are described by the activation energy losses $U_{\mathrm{act}}$ through

$$\begin{array}{c}\hfill i^{\prime \prime }=i_0^{\prime \prime }·\left[exp\left({\displaystyle \frac{\alpha ·z·F}{R_{\mathrm{m}}·T}}·U_{\mathrm{act}}\right)-exp\left({\displaystyle \frac{(1-\alpha )·z·F}{R_m·T}}·U_{\mathrm{act}}\right)\right]\end{array}$$

(2)

the ohmic losses $U_{\mathrm{ohm}}$ by

$$\begin{array}{cc}\hfill U_{\mathrm{ohm}}& =i^{\prime \prime }·A·\left(\sum R_{\mathrm{Elec}}+{\displaystyle \frac{{\delta }_{\mathrm{mem}}}{A·\kappa }}\right)\hfill \end{array}$$

(3)

and at last the mass transport losses $U_{\mathrm{diff}}$ using the following equation:

$$\begin{array}{cc}\hfill U_{\mathrm{diff}}& ={\displaystyle \frac{R·T}{z·F}}ln\left({\displaystyle \frac{p^0}{p^{*}}}\right)\hfill \end{array}$$

(4)

In the case of activation and mass transport losses, the calculations are made individually for the cathode and anode side. In addition, the effects depend on the exchange $i_0^{\prime \prime }$ and actual current density $i^{\prime \prime }$, the transfer coefficient of the reaction $\alpha$, the electric resistance $R_{\mathrm{Elec}}$, the membrane thickness ${\delta }_{\mathrm{mem}}$, the proton conductivity $\kappa$ and the pressure $p^{*}$ of $H_2$ and $O_2$ at the catalytic layer [45].

3. Methodology

The work follows the V-model for the development of control software, from the formulation of various requirements for the control model to system testing on a Rapid Control Prototyping (RCP) control unit [46]. In the course of the development process, the fuel cell control unit (FCCU) with the MPC is implemented in the existing overall software, and the control functionalities are validated in a MiL environment before deploying it on the RCP control unit.

The simulation environment consists of a control unit (CU) and a physical model, both implemented in MATLAB Simulink. The control units of the individual components are coordinated by a higher-level control unit, the hybrid control unit (HCU). The input variables in the MPC are the driving profile information of the heavy-duty truck consisting of the velocity, acceleration and gradient profile, a target SoC and the feed-back system state variables. The result of the optimization problem is transferred to the system model in the form of an optimized power trajectory of the fuel cell. For a functioning control loop, the environment of the embedded system model must be simulated. For this purpose, the system model consisting of the battery, vehicle and fuel cell is set up as a physical model within the MiL. An overview of the Optimal Control Problem (OCP) with the control constraints and the minimization objectives is shown in Figure 4.

3.1. Models of the FCEV

3.1.1. PEMFC Model

The MPC model used for the description of fuel cell dynamics is a 0-D fuel cell model based on the polarization curve of the fuel cell. The characteristic curve is defined as the operating voltage over the operating current, with the voltage losses being deducted, and is set up in the form of an equation. The Nernst voltage $U_{\mathrm{N}}$, also known as the open circuit voltage (OCV), is lower than the thermodynamically determined theoretical voltage value. Through a series of irreversible processes, the cell voltage is further reduced, resulting in the cell voltage $U_{\mathrm{A}}$.

$$\begin{array}{cc}\hfill U_{\mathrm{A}}& =U_{\mathrm{N}}-\Delta U\hfill \end{array}$$

(5)

$$\begin{array}{cc}\hfill \Delta U& =\Delta U_{\mathrm{Act}}+\Delta U_{\mathrm{Ohm}}+\Delta U_{\mathrm{Conc}}\hfill \end{array}$$

(6)

The voltage losses can be considered to consist of activation overvoltage in combination with electrochemical limitations at the electrode surfaces $\Delta U_{\mathrm{Act}}$, ohmic losses in both the ionic and electronic conductors $\Delta U_{\mathrm{Ohm}}$, and concentration overvoltages at the electrodes during the reaction process $\Delta U_{\mathrm{Conc}}$ [47]. By simplifying the fuel cell model using the polarization curve, the cell voltage $U_{\mathrm{A},\mathrm{BoL}}$ is approximated as a function of the stack current $I_{\mathrm{St}}$ resulting from the power demand $P_{\mathrm{fc}}$ of a single cell.

$$\begin{array}{cc}\hfill U_{\mathrm{A},\mathrm{BoL}}& =f\left(I_{\mathrm{St}}\right)\hfill \end{array}$$

(7)

$$\begin{array}{cc}\hfill I_{\mathrm{St}}& =f\left(P_{\mathrm{fc}}\right)\hfill \end{array}$$

(8)

The fuel cell efficiency ${\eta }_{\mathrm{Fc}}$ related to the lower calorific value is the ratio of the net electrical system output $P_{\mathrm{el}}$ to converted chemical energy and includes the fuel mass flow ${\dot{m}}_{{\mathrm{H}}_2}$ and the lower calorific value of hydrogen $H_{\mathrm{u},{\mathrm{H}}_2}$.

$$\begin{array}{c}\hfill {\eta }_{\mathrm{Fc}}={\displaystyle \frac{P_{\mathrm{el}}}{{\dot{m}}_{{\mathrm{H}}_2}·H_{\mathrm{u},{\mathrm{H}}_2}}}={\displaystyle \frac{U_{\mathrm{St}}·I_{\mathrm{St}}}{{\dot{m}}_{{\mathrm{H}}_2}·H_{\mathrm{u},{\mathrm{H}}_2}}}\end{array}$$

(9)

The stoichiometric hydrogen consumption ${\dot{m}}_{{\mathrm{H}}_2}$ can be calculated using Faraday’s law, where the molar mass of the hydrogen $M_{{\mathrm{H}}_2}$, the Faraday constant F and the electrochemical value z are used in the following equation:

$$\begin{array}{cc}\hfill {\dot{m}}_{{\mathrm{H}}_2}& ={\displaystyle \frac{M_{{\mathrm{H}}_2}·P_{\mathrm{el}}}{z·U_{\mathrm{A}}·F}}\hfill \end{array}$$

(10)

The power of the auxiliary loads, consisting of the components described in sec:SimModel, is determined as a function of the current fuel cell power. With an increase in the power demand on the fuel cell stack, there is a rise in the stack current, thereby determining the power of the BoP auxiliaries. The power output of the fuel cell system $P_{\mathrm{Fc},\mathrm{Sys}}$ is therefore the stack power $P_{\mathrm{Fc},\mathrm{St}}$ minus the power of the auxiliary loads $P_{\mathrm{BoP}}$.

$$\begin{array}{cc}\hfill P_{\mathrm{BoP}}& =f\left(I_{\mathrm{St}}\right)\hfill \end{array}$$

(11)

$$\begin{array}{cc}\hfill P_{\mathrm{Fc},\mathrm{Sys}}& =P_{\mathrm{Fc},\mathrm{St}}-P_{\mathrm{BoP}}\hfill \end{array}$$

(12)

In this work, the BoP power is modeled at system level as a static surrogate of the detailed plant model: the individual auxiliary consumers (air path, coolant circuit, etc.) are represented by map-based component models in the MiL plant, and their electrical power demands are summed up to obtain a total auxiliary power. For the EMS/MPC formulation, this total is reduced to a one-dimensional characteristic $P_{\mathrm{BoP}}=f\left(I_{\mathrm{St}}\right)$ (equivalently $f\left(P_{\mathrm{Fc},\mathrm{St}}\right)$), obtained from steady-state operating points of the detailed system model and implemented as a lookup table with interpolation. This choice keeps the OCP formulation simple and computationally efficient while still capturing the dominant trend that higher stack currents require higher auxiliary power.

In order to keep the degradation mechanisms of the fuel cell as low as possible, this work aims at a constant operation of the fuel cell with respect to the power output. Dynamic operation of the fuel cell leads to a high mechanical degradation of the fuel cell membrane and should therefore be avoided. In [48], a limitation of the voltage drop rate of the fuel cell is used to develop a strategy to protect against degradation. The aim of the strategy is to define an adaptive dynamic limit on the output power by formulating a relationship between the change in output power over time and the voltage drop rate. The drop rate ${\displaystyle \frac{\partial u_{\mathrm{decay}}}{\partial t}}$ is calculated based on [48] and [49] using the following equation for each node k on the prediction horizon:

$$\begin{array}{c}\hfill {\displaystyle \frac{\partial u_{\mathrm{decay}}}{\partial t}}=-{\displaystyle \frac{225\times \sigma [P_{\mathrm{Fc},\mathrm{St}}\left(k\right)+\dots +P_{\mathrm{Fc},\mathrm{St}}(k-4)]+10}{3600}}\end{array}$$

(13)

Experimental results in [48] show that the voltage drop rate under unit standard deviation $\sigma$ of five consecutive power steps is 225 ${\displaystyle \frac{\mu V}{h}}$ in dynamic operation and 10 ${\displaystyle \frac{\mu V}{h}}$ in steady-state operation. The same approach was used to determine the drop rate in [49]. This is an empirical model for the effects of cyclic operation on the fuel cell degradation. According to [49], the degradation mechanisms caused by cyclic operation are mechanical wear of the membrane due to fluctuations in temperature and potential, and chemical degradation due to fluctuations in humidity and detachment of platinum particles at the cathode. The voltage drop rate calculates the decrease of the possible open circuit voltage using the standard deviation for five consecutive fuel cell power demands $P_{\mathrm{Fc},\mathrm{St}}[k...k-4]$. The MPC calculates the state variables based on the optimal control trajectory for the predicted time horizon into the future. The change in fuel cell power over the past five time steps cannot be stored in the MPC calculation and therefore cannot be used to calculate a standard deviation as shown in [48] and [49]. Alternatively, the change in fuel cell power can be interpreted as the difference between two consecutive time steps and used as the basis for an approximation of the standard deviation. Using the calculation within the MPC, the standard deviation corresponds to the change in power from one time step to the next according to the following assumption:

$$\begin{array}{c}\hfill {\sigma }_{\mathrm{MPC}}={\displaystyle \frac{\frac{P_{\mathrm{Fc},\mathrm{St}}(k+1)-P_{\mathrm{Fc},\mathrm{St}}\left(k\right)}{P_{\mathrm{Fc},\max }-P_{\mathrm{Fc},\min }}}{2}}={\displaystyle \frac{\Delta P_{\mathrm{Fc},\mathrm{St}}}{2·(P_{\mathrm{Fc},\max }-P_{\mathrm{Fc},\min })}}\end{array}$$

(14)

The value of the approximate standard deviation is therefore defined between 0 and 1. At 0, the fuel cell is operating in steady state. If the value of the approximate standard deviation is 1, the fuel cell is operating in maximum dynamic mode and the cell voltage drop is equal to the maximum defined rate of decrease. This means that the voltage drop can be approximated and considered as an additional term within the cost function with respect to the degradation of the fuel cell system.

The actual cell voltage $U_{\mathrm{A},\mathrm{act}}$ can be calculated from the continuous decrease in cell voltage over time using the following equation:

$$\begin{array}{cc}\hfill U_{\mathrm{A},\mathrm{act}}& =U_{\mathrm{A},\mathrm{BoL}}+U_{\mathrm{decay}}\left(t\right)\hfill \end{array}$$

(15)

$$\begin{array}{cc}\hfill U_{\mathrm{decay}}\left(t\right)& ={\int }_{t_{\mathrm{BoL}}}^t{\displaystyle \frac{\partial u_{\mathrm{decay}}}{\partial t}}·dt\hfill \end{array}$$

(16)

The calculation of the decay rate can be used as a parameter to calculate the state of health (SoH) of the fuel cell. The SoH can be determined in the form of the following equation according to [49]:

$$\begin{array}{cc}\hfill \mathrm{SoH}& =1-{\displaystyle \frac{U_{\mathrm{A},\mathrm{BoL}}-U_{\mathrm{A}}}{U_{\mathrm{A},\mathrm{BoL}}-U_{\mathrm{A},\mathrm{EoL}}}}\hfill \end{array}$$

(17)

Here, $U_{\mathrm{A}}$ denotes the current cell voltage and $U_{\mathrm{A},\mathrm{EoL}}$ the cell voltage at end of life (EoL). Following [50], the EoL is reached as soon as the cell voltage at the rated operating point has dropped by 10 % with respect to its value at the beginning of life (BoL):

$$\begin{array}{cc}\hfill U_{\mathrm{A},\mathrm{EoL}}& =0.9·U_{\mathrm{A},\mathrm{BoL}}\hfill \end{array}$$

(18)

Based on , the SOH can be calculated as follows:

$$\begin{array}{cc}\hfill \mathrm{SoH}& =1-{\displaystyle \frac{U_{\mathrm{decay}}}{0.1}}\hfill \end{array}$$

(19)

3.1.2. Vehicle Model

The driving resistance and the resulting power requirement of the vehicle are determined through longitudinal dynamics with the velocity, acceleration and gradient profile of the driving cycles. The total driving resistance $F_t$ is composed of various partial resistances:

$$\begin{array}{cc}\hfill F_t\left(t\right)& =F_{\mathrm{air}}\left(t\right)+F_{\mathrm{r}}\left(t\right)+F_{\mathrm{g}}\left(t\right)+F_{\mathrm{a}}\left(t\right)\hfill \end{array}$$

(20)

$F_{\mathrm{air}}\left(t\right)$ is the aerodynamic resistance, which is based on the air friction on the vehicle surface and is given by

$$\begin{array}{cc}\hfill F_{air}\left(t\right)& ={\displaystyle \frac{1}{2}}·{\rho }_{\mathrm{air}}·A_{\mathrm{proj}}·c_{\mathrm{w}}·v{\left(t\right)}_{\mathrm{act}}^2\hfill \end{array}$$

(21)

where ${\rho }_{\mathrm{air}}$ is the air density, $A_{\mathrm{proj}}$ the projection area of the vehicle, $c_{\mathrm{w}}$ the aerodynamic drag coefficient and $v_{\mathrm{act}}$ the current velocity.

The rolling resistance $F_{\mathrm{r}}\left(t\right)$ is composed of the acceleration due to gravity g, the road slope angle $\alpha$, the vehicle mass $m_{\mathrm{v}}$ and a constant rolling resistance coefficient $f_{\mathrm{r}}$. It is calculated as follows:

$$\begin{array}{cc}\hfill F_{\mathrm{r}}\left(t\right)& =f_{\mathrm{r}}·m_{\mathrm{v}}·g·cos\left(\alpha \right)\hfill \end{array}$$

(22)

$F_{\mathrm{g}}\left(t\right)$ represents the grade force due to road slope, which is calculated as follows:

$$\begin{array}{cc}\hfill F_{\mathrm{g}}\left(t\right)& =m_{\mathrm{v}}·g·sin\left(\alpha \right)\hfill \end{array}$$

(23)

The resistance due to acceleration $F_{\mathrm{a}}\left(t\right)$ is calculated from the total mass $m_{\mathrm{v}}$ multiplied by the acceleration $a_{\mathrm{v}}$:

$$\begin{array}{cc}\hfill F_{\mathrm{a}}\left(t\right)& =m_{\mathrm{v}}·a_{\mathrm{v}}\left(t\right)\hfill \end{array}$$

(24)

Multiplying the total drag force on the wheel $F_{\mathrm{t}}$ with the current vehicle velocity $v_{\mathrm{act}}$ results in the required power $P_{\mathrm{t}}$ at the wheel:

$$\begin{array}{cc}\hfill P_{\mathrm{t}}\left(t\right)& =F_{\mathrm{t}}\left(t\right)·v_{\mathrm{act}}\hfill \end{array}$$

(25)

The base vehicle and dimensioning are based on [40]. The specifications of the truck are given in Table 1:

3.1.3. Electric Motor Model

The power of the electrical machine (EM) $P_{\mathrm{em}}$ can thus be determined via the efficiency chain, according to the following equation:

$$\begin{array}{c}\hfill P_{\mathrm{em}}\left(t\right)=\left\{\begin{array}{cc}{\displaystyle \frac{P_{\mathrm{t}}}{{\eta }_{\mathrm{dif}}·{\eta }_{\mathrm{gb}}·{\eta }_{\mathrm{inv}}}}+P_{{\mathrm{em}}_{\mathrm{loss}}},\hfill & P_{\mathrm{t}}\ge 0\hfill \\ P_{\mathrm{t}}·{\eta }_{\mathrm{dif}}·{\eta }_{\mathrm{gb}}·{\eta }_{\mathrm{inv}}+P_{{\mathrm{em}}_{\mathrm{loss}}},\hfill & P_{\mathrm{t}}<0\hfill \end{array}\right.\end{array}$$

(26)

where ${\eta }_{\mathrm{dif}}$ is the efficiency of the differential, ${\eta }_{\mathrm{gb}}$ is the efficiency of the gearbox and ${\eta }_{\mathrm{inv}}$ is the efficiency of the inverter.

The power loss of the electrical machine $P_{{\mathrm{em}}_{\mathrm{loss}}}$ is determined via a characteristic map as a function of the rotational speed $n_{\mathrm{em}}$ and the torque $T_{\mathrm{em}}$. Using a third-degree polynomial function the power loss can be approximated according to the following equation:

$$\begin{array}{c}\hfill \begin{array}{c}\hfill P_{{\mathrm{em}}_{\mathrm{loss}}}=a_0+a_1·n_{\mathrm{em}}+a_2·|T_{\mathrm{em}}|+a_3·n_{\mathrm{em}}·\left|T_{\mathrm{em}}\right|\\ \hfill +a_4·|T_{\mathrm{em}}{|}^2+a_5·n_{\mathrm{em}}·|T_{\mathrm{em}}{|}^2+a_6·{\left|T_{\mathrm{em}}\right|}^3\end{array}\end{array}$$

(27)

The maximum drive power $P_{{\mathrm{prop}}_{\max }}$, which is composed of the power of the energy sources and the power for the auxiliary consumers, can be defined as follows:

$$\begin{array}{c}\hfill P_{{\mathrm{prop}}_{\max }}={\eta }_{\mathrm{inv}}·\left(P_{{\mathrm{bat}}_{\max }}+\left(P_{{\mathrm{Fc}}_{\max }}-P_{\mathrm{BoP}}\right)·{\eta }_{{\mathrm{dcdc}}_{\mathrm{Fc}}}\right)\end{array}$$

(28)

3.1.4. Battery Model

For simplification, the focus is exclusively on the electrical battery model, as the thermal battery model is considerably more complex and challenging to represent due to the thermally dependent electrochemical processes [51,52]. In the electrical battery model, the voltage and the resulting current of a battery cell are calculated. The open circuit voltage of a battery cell $u_{\mathrm{OCV}}$ is determined based on the battery temperature $T_{\mathrm{bat}}$ and the $SoC$ characteristic field and described in the form of a first-degree interpolated function for the purpose of simplification.

$$\begin{array}{c}\hfill u_{\mathrm{OCV}}=3.262+0.892·SoC-9.044·{10}^{-18}·T_{\mathrm{bat}}\end{array}$$

(29)

The battery temperature $T_{\mathrm{bat}}$ is considered constant and is set to 25 °C In order to simplify the model, the thermal and transient temperature effects are neglected, since the battery model is now purely electric [53]. The only remaining state variable for the battery voltage is therefore the SoC. The high-voltage lithium-ion battery unit has a 186s3p configuration with a maximum capacity of 140 Ah and 98.84 kWh of energy. According to the specification, there are 186 cells connected in series $n_{\mathrm{bat},\mathrm{ser}}$ and 3 parallel groups of cells $n_{\mathrm{bat},\mathrm{par}}$. This results in the total battery voltage with open terminals $U_{\mathrm{OCV}}$ and the total capacity $C_{\mathrm{bat}}$ according to the following equations, with $C_{\mathrm{cell}}$ representing the capacity of a single battery cell:

$$\begin{array}{cc}\hfill U_{\mathrm{OCV}}& =u_{\mathrm{OCV}}·n_{\mathrm{bat},\mathrm{ser}}\hfill \end{array}$$

(30)

$$\begin{array}{cc}\hfill C_{\mathrm{bat}}& =C_{\mathrm{cell}}·n_{\mathrm{bat},\mathrm{par}}\hfill \end{array}$$

(31)

Following the subtraction of the internal resistance losses, the current battery voltage is determined. The loss voltage is contingent on the resulting current.

$$\begin{array}{c}\hfill U_{\mathrm{bat},\mathrm{act}}=U_{\mathrm{OCV}}-R_{\mathrm{int}}·I_{\mathrm{bat},\mathrm{act}}\end{array}$$

(32)

$$\begin{array}{c}\hfill R_{\mathrm{int}}=\left\{\begin{array}{cc}{R}_{\mathrm{int},\mathrm{dis}}·{\displaystyle \frac{n_{\mathrm{bat},\mathrm{ser}}}{n_{\mathrm{bat},\mathrm{par}}}},\hfill & I_{\mathrm{bat},\mathrm{act}}\ge 0\hfill \\ R_{\mathrm{int},\mathrm{chr}}·{\displaystyle \frac{n_{\mathrm{bat},\mathrm{ser}}}{n_{\mathrm{bat},\mathrm{par}}}},\hfill & I_{\mathrm{bat},\mathrm{act}}<0\hfill \end{array}\right.\end{array}$$

(33)

The battery model has been reduced to a static equivalent circuit diagram, as illustrated in Figure 5.

According to [54] and [13], the lithium-ion battery can be described by the internal electrical resistance model, using the following equation as a function of the voltage at the open terminals $U_{\mathrm{OCV}}$, the internal resistance $R_{\mathrm{int}}$, the power $P_{\mathrm{bat}}$ and the capacity of the battery $C_{\mathrm{bat}}$.

$$\begin{array}{c}\hfill {\displaystyle \frac{\partial \mathrm{SoC}}{\partial t}}=-{\displaystyle \frac{U_{\mathrm{OCV}}-\sqrt{U_{\mathrm{OCV}}^2-4·P_{\mathrm{bat}}·R_{\mathrm{int}}}}{2·R_{\mathrm{int}}·C_{\mathrm{bat}}}}\end{array}$$

(34)

From the relationship

$$\begin{array}{c}\hfill \mathrm{SoC}={\mathrm{SoC}}_0-{\displaystyle \frac{1}{C_{\mathrm{bat}}}}\int I_{\mathrm{bat},\mathrm{act}}dt\end{array}$$

(35)

with the initial state of charge ${\mathrm{SoC}}_0$ gives the actual current of the battery:

$$\begin{array}{c}\hfill I_{\mathrm{bat},\mathrm{act}}={\displaystyle \frac{U_{\mathrm{OCV}}-\sqrt{U_{\mathrm{OCV}}^2-4·P_{\mathrm{bat}}·R_{\mathrm{int}}}}{2·R_{\mathrm{int}}}}\end{array}$$

(36)

The current battery power $P_{\mathrm{bat},\mathrm{act}}$ is composed of the required power of the electric motor $P_{\mathrm{em}}$, the supply power of the components of the fuel cell $P_{\mathrm{BoP}}$ and the electrical power supplied to the fuel cell $P_{\mathrm{Fc}}$.

$$\begin{array}{c}\hfill P_{\mathrm{bat},\mathrm{act}}=P_{\mathrm{em}}-(P_{\mathrm{Fc}}-P_{\mathrm{BoP}})·{\eta }_{\mathrm{dcdc},\mathrm{Fc}}\end{array}$$

(37)

The maximum charge power $P_{\mathrm{bat},\min }$ and discharge power $P_{\mathrm{bat},\max }$ of the battery depends on the battery temperature $T_{\mathrm{bat}}$ and the SoC. The map-based value can be approximated by a polynomial function for a battery temperature of 25 °C:

$$\begin{array}{cc}\hfill P_{\mathrm{bat},\max }& =f\left(\mathrm{SoC}\right)\hfill \end{array}$$

(38)

$$\begin{array}{cc}\hfill P_{\mathrm{bat},\min }& =f\left(\mathrm{SoC}\right)\hfill \end{array}$$

(39)

The parameters of the components used in this study are given in Table 2.

3.2. MPC Energy Management Strategy

3.2.1. ACADOS Solver

acados is an open source software package that provides a collection of solvers for optimization problems. acados is based on a powerful linear algebra library Blasfeo with an interface to MATLAB and Python and is designed for fast optimization in embedded applications [38]. It uses the sequential quadratic programming (SQP) method to solve non-linear optimization problems (NOP). The Lagrange-Newton method for constrained problems with equation constraints is extended with inequality constraints and attempts to solve the optimization problem by iteratively solving quadratic approximations in the form of a linear Taylor approximation [55]. For this purpose, the Karush-Kuhn-Tucker (KKT) condition is used as an extension of the Lagrange method to solve the nonlinear optimization problem with inequality constraints and equation constraints.

3.2.2. Problem Formulation

To simplify the control problem, the change in fuel cell power is defined as the target variable. This has the advantages of avoiding excessive oscillations and of allowing the power gradient to be defined as a term within the cost function. The power of the battery is therefore equal to the difference between the power of the fuel cell and the power required by the powertrain and its consumers. The system input vector can be defined as follows:

$$\begin{array}{cc}\hfill u& =\left(\Delta P_{\mathrm{Fc}}\right)\hfill \end{array}$$

(40)

The state variables x consist of the SoC, the consumed hydrogen $m_{{\mathrm{H}}_2}$, the current fuel cell power output $P_{\mathrm{Fc}}$ and the cell voltage loss of the fuel cell $U_{\mathrm{decay}}$ of the dynamic system model:

$$\begin{array}{cc}\hfill x& =\left(\begin{array}{c}\mathrm{SoC}\\ m_{{\mathrm{H}}_2}\\ P_{\mathrm{Fc}}\\ U_{\mathrm{decay}}\end{array}\right)\hfill \end{array}$$

(41)

The dynamic behavior of the system can be expressed as follows:

$$\begin{array}{cc}\hfill \dot{x}& =\left(\begin{array}{c}-{\displaystyle \frac{I_{\mathrm{bat},\mathrm{act}}}{C_{\mathrm{bat}}·3600}}\\ {\displaystyle \frac{P_{\mathrm{Fc}}·M_{{\mathrm{H}}_2}}{2·F·U_{A,act}}}\\ \Delta P_{\mathrm{Fc}}\\ -{\displaystyle \frac{0.000225}{3600}}·\sigma -{\displaystyle \frac{0.000010}{3600}}\end{array}\right)\hfill \end{array}$$

(42)

System constraints help to take into account boundary conditions resulting from safety aspects or physical limitations of control and state variables. Hard constraints are defined for the state and control variables in the form of inequalities whose limits must not be exceeded. The following constraints are defined for the state and control variables:

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill {\mathrm{SoC}}_{\min }& \le \mathrm{SoC}\le {\mathrm{SoC}}_{\max }\hfill \\ \hfill 0& \le m_{H_2}\le m_{{\mathrm{H}}_2,\max }\hfill \\ \hfill P_{\mathrm{Fc},\min }& \le P_{\mathrm{Fc}}\le P_{\mathrm{Fc},\max }\hfill \\ \hfill U_{\mathrm{decay},\min }& \le U_{\mathrm{decay}}\le U_{\mathrm{decay},\max }\hfill \\ \hfill \Delta P_{\mathrm{Fc},\min }& \le \Delta P_{\mathrm{Fc}}\le \Delta P_{\mathrm{Fc},\max }\hfill \end{array}\end{array}$$

(43)

The SoC of the battery is limited to protect against overcharging or deep discharging and to protect against aging or performance degradation. The total hydrogen consumption is limited by the size of the tank. In addition to the physical limitation, hydrogen consumption is factored into the cost function and is therefore an optimization objective that should be minimized in the multivariable cost function with respect to minimum fuel consumption. Both the fuel cell and the battery have physical limits that must not be exceeded. A minimum power is defined as the lower limit for the fuel cell. This is to prevent the fuel cell from switching on and off and to prevent degradation. A soft limit is defined for the cell voltage of a fuel cell:

$$\begin{array}{cc}\hfill U_{\mathrm{A},\mathrm{act}}& \le U_{\mathrm{A},\max }+{ϵ}_1\hfill \end{array}$$

(44)

$$\begin{array}{cc}\hfill U_{\mathrm{A},\mathrm{act}}& \ge U_{\mathrm{A},\min }-{ϵ}_2\hfill \end{array}$$

(45)

The parameters ${ϵ}_1$, ${ϵ}_2$ allow the limits $U_{\mathrm{A},\min }$, $U_{\mathrm{A},\max }$ to be exceeded and are called slack variables. The slack variables appear within a term in the cost function and are penalized if the limit is exceeded. The purpose of limiting the cell voltage is to avoid excessively high and low voltages in the cell which increase cell degradation.

The control and state variables to be optimized are combined as cost terms within the cost function which is calculated as follows:

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill minJ\left(x\right(k),u(k\left)\right)& ={\displaystyle \sum _{k=0}^{N_{\mathrm{p}}-1}}(\beta {\left(k\right)}^{\mathrm{T}}·Q·\beta \left(k\right)·Y)\hfill \\ & +{\displaystyle \sum _{k=0}^{N_{\mathrm{p}}-1}}\left({\displaystyle \frac{1}{2}}·{ϵ}^{\mathrm{T}}·G·ϵ\right)\hfill \\ & +\left(A·{(\mathrm{SoC}\left(N_{\mathrm{p}}\right)-{\mathrm{SoC}}_{\mathrm{ref}}\left(N_{\mathrm{p}}\right))}^2·W\right)\hfill \end{array}\end{array}$$

(46)

with

$$\begin{array}{cc}\hfill \beta \left(k\right)& =\left(\begin{array}{c}{\dot{m}}_{{\mathrm{H}}_2}\left(k\right)\\ \sigma \left(k\right)\\ 1-{\eta }_{\mathrm{ec},\mathrm{Fc}}\left(k\right)\\ u\left(k\right)\end{array}\right)\hfill \end{array}$$

(47)

The cost function can be roughly divided into three terms. The first term is the vector $\beta \left(k\right)$ describing the objectives to be minimized within the cost function at each time step within the horizon H for the nodes $k=0\dots N_{\mathrm{p}-1}$. The second term is the penalty for exceeding the cell voltage limits of a fuel cell and the third term is the deviation of the SoC from a reference SoC value at the end of the horizon $k=N_{\mathrm{p}}$ as the terminal cost. Since only the SoC is given as a reference, this term describes the deviation of the SoC from the reference SoC. The matrices $Q,Y$ and the variables $A,W$ are responsible for factorizing and weighting of the respective terms within the cost function.

3.2.3. Drive Cycles

The VECTO Long Haul cycle is a driving cycle for determining the $C{O}_2$ emissions and the fuel consumption of heavy duty commercial vehicles. It is derived from the European Commission’s Vehicle Energy Consumption Calculation Tool (VECTO) [56]. The cycle represents a typical truck driving profile with constant speed sections and few start-stop operations over a distance of approximately 100 km with uphill and downhill gradients at a maximum speed of up to 87.5 km/h. Figure 6 shows the driving profile information consisting of acceleration, velocity and gradient profile over the cycle time.

Another cycle used to validate the MPC is a real world driving cycle from the city of Aachen through the Eifel region. The Eifel is an area in North Rhine-Westphalia, Germany, which provides an interesting cycle for comparing the MPC with the DP due to its dynamic elevation and velocity profile. The acceleration, velocity and gradient profile over the cycle time is shown in Figure 7.

4. Results

The fuel cell degradation reported in this section is calculated by an ECSA-aging model integrated into the stack model described in sec:SimModel, following established literature on electrochemically active surface area (ECSA) aging [57,58,59,60]. Within the MPC cost function, a cyclic-operation-induced voltage decay is represented by the state $U_{\mathrm{decay}}$, obtained by integrating the empirically parameterized decay rate (see 15) with the dynamic-stress proxy $\sigma$ (see 14), so that the effect of power fluctuations on the voltage loss is directly penalized during optimization. The resulting end-of-cycle degradation is evaluated through two complementary indicators based on the aging model from Thiele et al. [57]: the relative membrane thickness change ($\Delta {\delta }_{\mathrm{mem}}$) and the relative electrochemically active surface area loss ($\Delta A_{\mathrm{ECSA}}$), both normalized to the DP strategy in Figure 12. While the absolute changes over a single drive cycle are small, the comparison highlights a consistent trend: the highly dynamic RB strategy leads to the largest degradation in both indicators and on both cycles, whereas the MPC attains the lowest membrane thickness loss and stays close to the globally optimal DP in terms of ECSA loss on the VECTO cycle, which is consistent with the objective of smoothing the fuel cell power trajectory.

4.1. Computing Time Analysis

In order to evaluate how the length of the prediction horizon H affects the computation time and the optimality of the results, a variation of the horizon lengths is performed and analyzed. Previously, the change in execution time $t_{\mathrm{CPU}}$ was investigated by varying the maximum iteration steps, with the result that the optimum solution was achieved within the defined maximum SQP iteration steps in each computation step, so this consideration is ignored. For this purpose, a cycle similar to the VECTO long haul is evaluated within the MiL for a horizon range of 50–1000 s. The considered execution time $t_{\mathrm{CPU}}$ is the total execution time of the S-function of the OCP in Simulink per sample time1. The maximum computation value occurring within the entire cycle is used as the reference value for the analysis, as this value is more meaningful in terms of execution on an ECU. The choice of prediction horizon has a direct influence on the execution time of the MPC. As the horizon increases, the optimum performance trajectory for each grid point is calculated for that period into the future. To quantify this influence, the cycle was run with the horizon lengths:

$$\begin{array}{c}\hfill H\in \{50\mathrm{s},100\mathrm{s},200\mathrm{s},300\mathrm{s},400\mathrm{s},500\mathrm{s},600\mathrm{s},800\mathrm{s},1000\mathrm{s}\}\end{array}$$

(48)

Figure 8 shows the progression of the target variables as a function of the horizon length of the MPC. It can be seen that the maximum execution time increases as the horizon length increases. The maximum execution time at 500 s and 600 s falls out of the pattern and is less than the execution time at 400 s. From a horizon length of 800 s, the maximum permissible execution time of 300 ms is exceeded. If executed on a controller, the MPC would run into task overruns, which must be avoided. Since the average execution time increases only slightly with increasing horizon length, but the cell voltage loss decreases significantly with increasing horizon length. At horizon length below 200 s the gradient of the cell voltage loss is significantly higher, thus with regards to the cell voltage loss a horizon length from at least 200 s should be chosen. It can be seen that a horizon length of 400 s reflects a good compromise between the lowest possible maximum execution time, low fuel consumption, low deviation from the target SoC and cell voltage loss. As a result of this parameter analysis, a prediction horizon of 400 s is selected for further investigation.

4.2. MPC vs. DP

4.2.1. VECTO Long Haul Cycle

In the following, the model predictive control approach is compared to a DP with respect to its optimality. The DP is designed to compute an optimal SoC trajectory, taking into account the minimum fuel consumption and the rate of change of fuel cell system power. The same boundary conditions and performance constraints are applied to both the MPC and the DP in order to compare the results. Due to the deviation of the mathematical system model in the MPC and DP from the simulation model in the MiL environment, a balanced state of charge at the end of the driving cycle cannot be guaranteed. In order to be able to guarantee the comparability of the consumption results nevertheless, the charge deviations in fuel consumption must be taken into account. A correction is therefore applied in both directions, so that all strategies are evaluated at an identical, balanced charge level at the end of the cycle. In the case of a positive charge deviation (${\mathrm{SoC}}_{\mathrm{end}}<{\mathrm{SoC}}_{\mathrm{start}}$), the energy missing in the battery is converted into an additional fuel mass to be applied, taking into account the overall efficiency, so that no artificial fuel consumption advantage is created. Conversely, in the case of a negative charge deviation (${\mathrm{SoC}}_{\mathrm{end}}>{\mathrm{SoC}}_{\mathrm{start}}$), which can occur, for example, due to extended recuperation during downhill driving, the surplus energy remaining in the battery is credited as an equivalent fuel mass and deducted from the consumption. Both cases are covered by the same correction rule formulated in [61] according to equation

$$\begin{array}{c}\hfill m_{{\mathrm{H}}_2,\mathrm{cor}}=m_{{\mathrm{H}}_2,\mathrm{end}}+{\displaystyle \frac{({\mathrm{SoC}}_{\mathrm{start}}-{\mathrm{SoC}}_{\mathrm{end}})·(-C_{\mathrm{batt}})}{{\eta }_{\mathrm{Fc}}·{\eta }_{\mathrm{inv}}·{\eta }_{\mathrm{bat}}}}·{\displaystyle \frac{1}{H_{\mathrm{u},{\mathrm{H}}_2}}}\end{array}$$

(49)

which adjusts the fuel mass $m_{{\mathrm{H}}_2,\mathrm{end}}$ determined from the simulation as a function of the charge difference. For the fuel cell system efficiency ${\eta }_{\mathrm{Fc}}$ an average efficiency from the drive cycle is used. The overall efficiency is the efficiency of the electrical path and includes the fuel cell, the inverter and the battery. The results of the corrected hydrogen consumption between the DP, MPC and RB strategy of the two drive cycles can be seen in Figure 11 and Figure 11.

The fuel cell system power and state of charge (SoC) trajectories of DP and MPC in the MiL environment using the gradient and speed profile for the VECTO Long Haul cycle are compared in Figure 9. The DP controller maintains a nearly constant fuel cell power of about 95 kW, which is close to one-third of the fuel cell’s maximum power and therefore in the range of the best efficiency point for fuel cell system operation. The power trajectory of the MPC, on the other hand, is considerably more dynamic and varies between a minimum of about 25 kW and a maximum of about 217 kW. These differences in power trajectories also affect the SoC trajectories. The DP allows a wider range of SoC to be covered compared to the MPC, due to its knowledge of the entire driving profile.

Within the prediction horizon $H=400\mathrm{s}$, the MPC searches for the optimal fuel cell power from which the SoC trajectory is derived. It tries to achieve the target SoC at the end of the horizon which is reflected in the maximum deviation of the SoC within the prediction horizon. The SoC of the DP reaches its minimum at $t=1770\mathrm{s}$ and drops to a value of about 41 %. The DP thus uses 34 % of the allowed SoC operating range (i.e., the permissible SoC band between ${\mathrm{SoC}}_{\min }$ and ${\mathrm{SoC}}_{\max }$). The SoC of the MPC reaches its minimum of 51 % at $t=1692\mathrm{s}$ and its global maximum of 66.8 % at $t=2094\mathrm{s}$, thus using 26 % of the allowed SoC range. From the results obtained, a SoC deviation from the desired target SoC can be determined for both the MPC and the DP. The MPC misses the target SoC by $-1.06\%$ and has a state of charge of 59.4 % at the end of the cycle, while the DP ends at a state of charge of 60.6 %, i.e. 0.93 % above the target SoC. After correcting the fuel consumption for the state of charge deviation at the end of the cycle, the fuel consumption for the MPC is 7.18 kg. With a fuel consumption of 6.97 kg, the globally optimal DP control therefore consumes 2.9 % less hydrogen than the MPC control.

4.2.2. Aachen-Eifel Real Driving Cycle

Figure 10 shows a comparison between the FC system power trajectories and the SoC trajectories of DP and MPC in a MiL environment using the gradient and speed profile of the Aachen-Eifel real driving cycle over time. The main difference between global and local optimal control is evident in the SoC power trajectories. The DP control, which benefits from global information about the driving profile, mostly performs in the range of optimal efficiency and minimum fuel consumption, similar to the VECTO Long Haul cycle. However, the DP control mostly operates in a narrow band between about 55 kW and 65 kW, with isolated peaks up to 95 kW due to the influence of the gradient and speed profile.

Unlike the MPC, the DP detects the high-speed section at the end of the cycle and charges the battery more initially, reaching a maximum SoC of 75.4 % at $t=1755\mathrm{s}$ before drawing it down to a minimum of 29.1 % at $t=4006\mathrm{s}$. As a result, the DP regulates a more constant fuel cell power output throughout the cycle to maintain low hydrogen consumption and high efficiency. On the other hand, the MPC reaches a global maximum SoC of 73.1 % at $t=4385\mathrm{s}$. Due to its limited knowledge of the driving cycle information according to the defined finite time horizon, the power output of the MPC differs from the power output of the DP, resulting in a more dynamic trajectory and a maximum power of 244 kW at $t=3865\mathrm{s}$. The DP uses a substantially larger part of the allowed SoC range (77 %) compared to the MPC (43 %). The final SoC of the DP control is 59.9 %, which is 0.2 % below the target SoC and the final SoC of the MPC control is 63.6 %, which is 6.1 % above the target SoC. With respect to the total hydrogen consumption, the advantage of the global optimality of the DP is clearly visible. After correcting the charge deviation for the MPC and DP control to the target SoC, the total hydrogen consumption of the MPC is 7.65 kg, which is 10.7 % higher than the DP’s consumption of 6.91 kg. Although the MPC does not reach the global optimum of the DP, it provides the power split in real time. In the next step the MPC is compared to a rule-based strategy.

4.3. MPC vs. Rule-Based Strategy

4.3.1. VECTO Long Haul Cycle

In this section, the MPC is compared to a rule-based (RB) strategy with respect to fuel consumption and end-of-cycle fuel cell stack degradation for the VECTO Long Haul drive cycle. The RB strategy corrects the fuel cell power demand based on the deviation of actual and reference SoC, taking into account the minimum and maximum defined fuel cell power limitations. The results of the desired fuel cell power and the SoC trajectory can be seen in Figure 9. Due to the penalization of the SoC deviation in the RB strategy, the power trajectory of the fuel cell is very dynamic with high power peak demands to minimize the deviation at all times. The power output varies between a minimum of about 25 kW and a maximum of 254.6 kW at $t=1739\mathrm{s}$. This can also be seen in the trajectory of the SoC with a global minimum SoC of 54.7 % at $t=1692\mathrm{s}$ and a maximum SoC of 66.3 % due to recuperation at about $t=2093\mathrm{s}$. Due to the non-predictive characteristic of the RB strategy and the braking process at the end of the cycle, a final SoC above the target SoC of 60 % is achieved. The final SoC of the RB strategy is 63.9 %, which is 6.4 % above the target SoC. After correcting the fuel consumption for the state of charge deviation at the end of the cycle, the total hydrogen consumption of the RB strategy is 7.67 kg, which is 6.9 % higher than the MPC’s consumption of 7.18 kg.

The end-of-cycle fuel cell degradation is evaluated via the relative membrane thickness change ($\Delta {\delta }_{\mathrm{mem}}$) and the relative electrochemically active surface area loss ($\Delta A_{\mathrm{ECSA}}$), based on the aging model from Thiele et al. [57] and normalized to the DP strategy in Figure 12. For the VECTO cycle, the dynamic RB strategy exhibits the highest degradation, with a membrane thickness loss 4.2 % and an ECSA loss 4.2 % above the DP level. The MPC, by contrast, remains close to the globally optimal DP, with a membrane thickness loss 2.1 % below and an ECSA loss only 1.8 % above the DP level, reflecting its smoother fuel cell power trajectory. The corresponding comparison for the RDE cycle is shown in Figure 12.

4.3.2. Aachen-Eifel Real Driving Cycle

In this section, the MPC is compared to a rule-based (RB) strategy with respect to fuel consumption and end-of-cycle fuel cell stack degradation for the Aachen-Eifel real driving cycle. The results of the desired fuel cell power and the SoC trajectory can be seen in Figure 10. Due to the more dynamic elevation and velocity profile, this is also reflected in the power output of the RB strategy with more fuel cell power peaks to minimize the SoC deviation at all times. This can also be seen in the trajectory of the SoC with a maximum charge deviation down to an SoC of 48.6 %. The SoC reaches a maximum value of about 74.4 % at $t=5338\mathrm{s}$. The final SoC of the RB strategy is 68.1 %, which is 13.5 % above the target SoC. After correcting the charge deviation for the RB strategy to the target SoC, the total hydrogen consumption is 8.75 kg, which is 14.3 % higher than the MPC’s consumption of 7.65 kg. The same trend is observed for the Aachen-Eifel real driving cycle in Figure 12. The RB strategy again shows the highest degradation, with a membrane thickness loss 3.5 % and an ECSA loss 4.1 % above the DP level. The MPC attains the lowest membrane thickness loss of all three strategies, 9.1 % below the DP level, but reaches an ECSA loss of 4.1 % above the DP level. Owing to its smoother power trajectory, the MPC therefore avoids the additional membrane degradation caused by the highly dynamic operation of the RB strategy, whereas in terms of ECSA loss the MPC and the RB strategy are practically on the same level on this cycle.

As shown in Figure 12, the DP achieves the lowest ECSA loss of the three strategies, whereas the MPC exhibits the lowest membrane thickness loss. This can be explained by the fact that the DP predominantly operates at a highly efficient, low-current operating point and therefore at a comparatively high cell voltage, whereas the MPC, due to its limited prediction horizon, more frequently commands higher load points, i.e., lower cell voltages, in order to reach the target SoC at the end of each horizon. Since the formation rate of hydrogen peroxide and the resulting chemical, voltage-dependent membrane degradation (e.g., via the associated fluoride release rate) is known to increase with cell voltage [62], higher operating voltages tend to accelerate membrane degradation. This voltage-dependent chemical degradation pathway is not implemented in the fuel cell model used within the hybrid energy management system, in order to limit modeling complexity (see [63] for a broader discussion of prescriptive lifetime management approaches for PEM fuel cell systems). This explains why the DP, despite its higher operating efficiency, does not also achieve the lowest membrane thickness loss. It should be noted that, as a first step, the aging models used in this work have only been checked for plausibility and have not yet been validated against real measurement data.

To assess the transferability of these findings to a different driving profile, all three hybrid strategies were additionally evaluated on the more dynamic Aachen-Eifel real driving cycle (cf. ). The qualitative trends observed for the VECTO Long Haul cycle are confirmed: the DP and the MPC end the cycle close to the target SoC, whereas the non-predictive RB strategy overshoots it by 13.5 %, and the ranking of the aging indicators follows the same pattern as for the VECTO cycle. However, the results also reveal that the advantage of the DP over the MPC becomes more pronounced under the more dynamic load profile of the Aachen-Eifel cycle, whereas the MPC achieves a closer approximation of the DP’s global optimum on the comparatively steadier VECTO Long Haul cycle.

Figure 11. Comparison of the corrected hydrogen consumption between the DP, MPC and RB strategy for the (fig:CompHydConsumption) VECTO Long Haul cycle and the (fig:CompHydConsumptionRural) Aachen-Eifel real driving (rural road) cycle.

Figure 11. Comparison of the corrected hydrogen consumption between the DP, MPC and RB strategy for the (fig:CompHydConsumption) VECTO Long Haul cycle and the (fig:CompHydConsumptionRural) Aachen-Eifel real driving (rural road) cycle.

Regarding applicability to a degraded stack, the methodology remains applicable as long as the MPC prediction model is updated to the actual stack condition. In practice, this can be achieved by (i) initializing the degradation-related state $U_{\mathrm{decay}}$ (and/or SoH) to the current value and (ii) updating the polarization-curve-based voltage model (and, if necessary, the efficiency/BoP maps) to reflect the degraded performance. The OCP structure, constraints, and real-time implementation remain unchanged; only the model parameters and initial conditions need to be adapted.

Figure 12. Comparison of the relative membrane thickness ($\Delta {\delta }_{\mathrm{mem}}$) and the relative electrochemically active surface area ($\Delta A_{\mathrm{ECSA}}$) at the end of the cycle between the DP, MPC and RB strategy for the (fig:CompAging) VECTO Long Haul cycle and the (fig:CompAgingRural) Aachen-Eifel real driving (rural road) cycle.

Figure 12. Comparison of the relative membrane thickness ($\Delta {\delta }_{\mathrm{mem}}$) and the relative electrochemically active surface area ($\Delta A_{\mathrm{ECSA}}$) at the end of the cycle between the DP, MPC and RB strategy for the (fig:CompAging) VECTO Long Haul cycle and the (fig:CompAgingRural) Aachen-Eifel real driving (rural road) cycle.

5. Conclusions

In this work, a model predictive control approach was developed for the calculation of an optimal fuel cell power trajectory for a heavy duty fuel cell truck. To compute the power trajectory, the optimization problem was formulated using the open source software package acados in MATLAB Simulink. Based on the cost function defined in this work, the optimal control trajectory is computed and transferred to the system model as a target fuel cell power curve. By taking into account the power gradient in the cost function, the MPC was able to achieve constant operation in terms of avoiding large power fluctuations, which has a positive effect on the lifetime of the fuel cell system. A good trade-off between computational time and optimality was found based on a parameter analysis of important target variables of the MPC to evaluate the influence of the length of the prediction horizon.

The developed optimization algorithm was successfully implemented in an existing model-in-the-loop environment and the functionality of the optimization algorithm was validated and tested using different driving cycles. In addition, a comparison of the developed MPC with a global optimal DP was performed with respect to hydrogen consumption and with a RB strategy with respect to aging of the fuel cell. The result of the comparison shows for typical truck driving profiles for the VECTO Long Haul a consumption advantage of the DP over the MPC of 2.9 %. For a more dynamic cycle that does not have a typical heavy duty truck driving profile, the MPC approach consumed 10.7 % more hydrogen than the DP for the Aachen-Eifel real driving cycle. The SoH of the fuel cell was derived by calculating the cell voltage drop based on the standard deviation of five consecutive fuel cell outputs over time as well as via the membrane thickness reduction. However, since the SoH depends on a large number of factors, the cell voltage loss calculated in this work based on empirically determined values was used for a simplified representation of the SoH. Due to the limited information about the fuel cell dynamics from the polarization curve, the simplified SoH is only intended to illustrate the degradation effects of power variations and does not represent the magnitude of various other effects. Within this work, the MPC reduced the end-of-cycle fuel cell degradation relative to the rule-based strategy. On the VECTO Long Haul cycle, the membrane thickness loss and the ECSA loss of the MPC were 6.1 % and 2.3 % lower than those of the RB strategy. On the dynamic RDE cycle the membrane thickness loss of the MPC was 12.2 % lower, whereas the ECSA loss of both strategies was practically identical.

With respect to different powertrain designs, the proposed control strategy is flexible because the OCP formulation is based on generic power and energy balances, component maps, and constraints. A change in component sizing (e.g., stack rated power, battery capacity, DC/DC limits, or auxiliary sizing) primarily affects the steady-state maps (fuel cell efficiency, BoP power) and the admissible operating region (power and SoC constraints), while the MPC structure and real-time implementation remain unchanged. This also makes the approach suitable as a control layer in co-optimization environments: additional objectives (e.g., thermal limits, operating cost, or mission-level constraints) and updated component models can be integrated by extending the stage cost and constraints without changing the underlying solver workflow. Finally, the general MPC/OCP framework is not restricted to heavy-duty road transport and can be adapted to other application domains (e.g., aviation) by replacing the vehicle/mission model, updating constraints and safety requirements, and using domain-specific fuel cell system models.

In the simulation environment, the specified cycle is assumed to be a real driving scenario in which the truck follows the speed profile. In real driving scenarios, this profile cannot be followed exactly due to various environmental influences. To capture the real driving situation within the MPC, telemetry data of the traffic and a predictive model, e.g. in the form of a Markov chain model, would be required to predict the new driving profile information and make it available to the MPC. By adding telemetry and a predictive model, MPC control could be used as a real-time application in real driving scenarios. In addition, the SoC trajectory generated by the DP could be used as a reference SoC for MPC control.