1. Introduction
Proton Exchange Membrane (PEM) fuel cells were first developed by General Electric in the United States in the 1960s and used in a spacecraft developed by the National Aeronautics and Space Administration (NASA) in the Gemini program [
1]. They are a type of acid fuel cell in which a solid or quasi-solid “membrane” material is used instead of the proton-conducting liquid electrolyte necessary in earlier acid cells. They can be used in cars, scooters, bicycles, boats, underwater vessels [
2,
3] and also aircrafts [
4]. These fuel cells have several advantages. One of them is low operation temperature which gives a fast start-up. Furthermore, a solid electrolyte is used so no electrolyte leakage can happen. PEM fuel cell design is also very simple and compact. Finally, they are very reliable in operation. On the other hand, there are minor drawbacks, such as the formation of a significant amount of water that may not have time to evaporate from the cell and thus damage it and the platinum catalyst’s high susceptibility to carbon monoxide. However, all of the positive features of the PEM fuel cells make them great alternatives to conventional power generation methods that rely on fossil fuels. Technological details related to different aspects of PEM fuel cells can be found in numerous publications, e.g., [
5,
6,
7].
Since PEM fuel cells are nonlinear processes, their efficient control is quite difficult. Different control strategies are possible for optimising the power split between the battery and PEM fuel cell to maximise system efficiency and reduce fuel consumption [
8]. Simple controllers, including a Sliding-Mode Controller (SMC) and a fuzzy logic controller, can supply a convenient voltage [
9]. A more advanced approach is to use a fuzzy Proportional-Integral-Derivative (PID) controller [
10] for temperature control of the fuel cell cooling system or an adaptive SMC algorithm which is integrated with the radial basis function neural network to control the current or voltage of the fuel cell power supply [
11]. Recent works show that PEM output power quality can be improved using classical and fractional-order PID controllers [
12]. An example of a fuzzy control method to control the temperature and humidity of the stack has been presented in [
13]. It is also possible to use a fuzzy adaptive PI controller with an improved advanced genetic algorithm [
14]. Hierarchical techniques can also be applied to develop a comprehensive energy management strategy consisting of high-level control, which is based on fuzzy logic technique as well as adaptive control loop and low-level control that generates a pulse-width-modulation signal so that appropriate power distribution is achieved and the operating voltage of the powertrain is stabilised [
15]. In recent years, applications of MPC algorithms [
16,
17,
18] to PEM fuel cells have been reported. A linear model can be used for MPC prediction [
19]. However, typically, a nonlinear model is used, meaning that a nonlinear optimisation task must be solved at each sampling instant on-line [
20,
21,
22,
23]. A practical method leading to simplification of the nonlinear task is successive on-line linearisation which leads to a quadratic optimisation problem [
24,
25].
In all MPC algorithms cited above, the MPC cost-function used considers the sum of absolute values of the predicted control errors (the L
norm). However, in a few publications concerned with other applications, it is suggested that better control is possible when the sum of absolute predicted control errors is utilised as the MPC cost function (the L
norm) [
26,
27,
28,
29,
30]. This paper aims to present a computationally simple MPC-L
control strategy and its application to the PEM process. This work extends previous research concerned with computationally efficient nonlinear MPC-L
in which a neural approximator of the absolute value function is utilised [
30]. Unfortunately, such a control strategy’s effectiveness heavily depends on the approximator’s accuracy. A more straightforward approach is recommended in this work. The neural network is replaced by two analytical approximators that are defined explicitly. The predicted process trajectory is linearised on-line to find a quadratic optimisation task. All implementation details of the discussed algorithm are detailed. Furthermore, the algorithm is thoroughly compared with the classical MPC-L
approach. A multi-criteria control quality assessment is performed as the MPC-L
and MPC-L
algorithms are compared using four control quality indicators. It is shown that the presented MPC-L
scheme gives better results for the PEM process.
The paper is organised as follows.
Section 2 presents the continuous-time PEM fuel system model used as the simulated process and a description of the Wiener model that is utilised in all MPC algorithms is next presented.
Section 3 recalls the general MPC optimisation task with both L
and L
norms. The main contribution is detailed in
Section 4 where a computationally efficient MPC-L
algorithm with two analytical approximators is described. The simulation results are presented and discussed in
Section 5. Finally, the research findings are summarised in
Section 6.
2. Proton Exchange Membrane Fuel Cell
Fuel cells are devices that allow electricity production, typically through the chemical conversion of hydrogen and oxygen. Fuel (hydrogen) is supplied continuously to the anode; oxygen is supplied to the cathode. An electrolyte differs depending on the fuel cell type. Energy and heat are generated during the electrochemical reaction; water is the only byproduct. In cars, the entire process begins with the supply of hydrogen from a high-pressure tank to the cell. In parallel, compressed air is also supplied. The result of the reaction in the cell is a generated current, which is converted into alternating current and supplied to the electric motor responsible for traction.
In this work, a PEM fuel cell is considered, which has the polymer membrane as the electrolyte that conducts protons. A catalyst covers the anode, facilitating the splitting of hydrogen atoms into protons and electrons. Ions pass through a polymer membrane conducting protons, which is an isolator for the electrons. Electrons flow through an external circuit, creating an external current in the cell. These fuel cells have been incredibly popular lately as alternative power sources. It is beneficial to use them due to their low operating temperature, great efficiency and minimal emission of harmful greenhouse gases.
This section describes a fundamental model of the PEM process. It is used in simulations only as the process. Next, a neural Wiener model is detailed. It is used in MPC algorithms only.
2.1. Description of the PEM Fuel Cell System
The PEM fuel cell model considered in this paper was first introduced in [
31] and further discussed in [
32,
33,
34]. The manipulated variable is the methane flow rate denoted as
q (
) and the output is the stack output voltage denoted as
V (V). The system is influenced by one measured disturbance which is the external current load, denoted as
I (A). The partial pressures of hydrogen, oxygen, and water are represented by the variables
,
and
, respectively (atm). The input hydrogen flow, hydrogen reacted flow, and oxygen input flow are denoted by
,
and
, respectively (
).
2.2. PEM Fuel Cell Continuous-Time Fundamental Model
The fundamental continuous-time model of the PEM system is described by the following continuous-time equations
Equation (
1) represents the hydrogen flow which can be obtained from the reformer, where
q is the methane flow rate while
,
and
are constants. The pressure of hydrogen is calculated from Equation (2), where
and
denote the valve molar constant for hydrogen and the response time of hydrogen flow, respectively. The current pressure of oxygen is found using Equation (3), where
and
denote the valve molar constant for oxygen and the response time of oxygen flow, respectively. Equations (4) and (5) describe the hydrogen flow that reacts and the pressure of water, respectively, where
is a constant. The stack voltage is found from Equation (6), where
E is calculated from the Nernst’s formula in Equation (7). Additionally, the following symbols are used:
is the number of cells in series in the stack,
is the ideal standard potential,
is the universal gas constant,
is the absolute temperature and
is the Faraday’s constant. Finally, activation losses and ohmic losses are described by Equations (8) and (9), respectively, where
is the internal resistance while
B and
C are activation voltage constants. The values of model parameters are given in
Table 1.
The values of process input and disturbance signals are constrained as follows
The values of process variables for the initial operating point are the following: , , , , and .
2.3. Wiener Model of PEM Fuel Cell
The discussed PEM fuel cell’s continuous-time fundamental model has linear dynamic transfer functions in Equations (2), (3) and (5), but the stack voltage is defined by the nonlinear steady-state Nernst’s equation (7). This means that the outputs of the linear dynamic part of the model are inputs of the nonlinear steady-state one. As a result, it is evident that a Wiener structure [
35] as an empirical model of the PEM fuel cell under consideration can be used. For the Wiener model, the process variables
q,
I, and
V are scaled as follows
where
,
and
denote values of the signals at the initial operating point.
Identification of a few variants of the Wiener model and selection of the best one is thoroughly discussed in [
24,
35]. The structure of the chosen model is shown graphically in
Figure 1. It has three linear dynamic blocks and one nonlinear steady-state block. The outputs of each linear block, at every sampling instant
k, are found with the use of the following formulas
The integers
for
,
for
,
and
define the order of the model dynamics. The constant parameters of the linear dynamic blocks are denoted by the real numbers
(
,
),
(
,
,
) and
(
). It is to be noted that the steady-state block of the model does not use just the measured outputs of the linear blocks,
,
,
, but also the measured disturbance signal,
h. The output of the model,
y, is calculated from the nonlinear steady-state block by using the following general equation
A neural network with four inputs, one hidden layer with
K units and one output is used as the nonlinear part of the model. The model output is
Weights of the network are denoted by , , and , , for the first and the second layers, respectively.
3. Problem Formulation
Although this work is concerned with a specific process (the PEM fuel cell), let us use universal notation in derivation of the MPC algorithm. The manipulated variable (the input of the process) is denoted as
u and the controlled variable (the output of the process) is denoted as
y. For scaling, Equation (
12) is used. Then, the classical rudimentary MPC optimisation task can be presented in the following form
The task is solved at each sampling instant and the solution is the vector of decision variables, of length equal to the control horizon denoted as
, containing the increments of the manipulated variable
The solution is found with respect to the constraints imposed on the magnitude of the input, defined by
and
, and the constraints imposed on the increments of that variable, defined by
and
. Only the first element of the decision variables vector is applied to the process, i.e., PEM fuel cell system. In general form, the minimised cost-function
can be written in the task as a sum of two parts as follows
The first part of the MPC cost-function is of the following form
The predicted control error for the future sampling instant
computed at the current time step
k is denoted by
. In the cost-function, they are taken into account over the prediction horizon
N as arguments of a function
. The predicted control errors are calculated as follows
where
is the set-point for the future sampling instant
available at the current time step
k and
is the predicted output of the process the sampling instant
found from a process model at the time step
k. The second part of the MPC cost-function is of the following form
All changes of the manipulated variable (over the control horizon) are minimised when this part is included into the cost-function
. In this part, there is also a weighting parameter
(the larger its value, the larger the penalty term and the slower the control). The function
that weighs the predicted control errors (Equation (
21)) can be a squared function
or it can be an absolute value function
Hence, if the square function (
24) is used, we obtain the classical MPC cost-function in which squared values of the predicted control errors are considered
It is called the L
cost-function and the resulting MPC algorithm is named MPC-L
. In this work, the L
cost-function is studied, which uses the part given by the Equation (
25). As a result, we obtain the following L
cost-function in which absolute values of the predicted control errors are considered
Let us note that both MPC cost-functions have the same second part defined by Equation (
23).
4. Computationally Efficient Nonlinear MPC Using the L Cost-Function
The nonlinear MPC-L
optimisation task is formulated with the use of the general form of the optimisation task in (
18) and the MPC-L
cost function (
27). It takes the following form
It is important to note that the PEM fuel cell system model is nonlinear, so the predictions of the process output, , are nonlinear functions of the calculated moves as well. In such a case, we have to distinguish two problems to solve.
Problem 1: The absolute value function used in the cost-function is not differentiable. As a result, the MPC-L cost-function is also not differentiable. It makes it impossible to use the classical gradient-based optimisation method. For the L norm, this problem does not exist because the MPC-L cost-function is quadratic and differentiable.
Solution to problem 1: The first part of the non-differentiable cost function (
27) is replaced by its differentiable representation. This work uses two analytical approximations of the absolute value function. It is also possible to solve this problem using neural approximators [
30].
Problem 2: Assuming that a nonlinear model is used in MPC, predictions
are nonlinear functions of the calculated decision vector (
19), which means that the MPC-L
or MPC-L
cost function is nonlinear. As a result, a nonlinear optimisation task is obtained and it must be repeatedly solved at each sampling instant.
Solution to problem 2: The
cost function after applying an analytical approximator is still nonlinear in terms of the computed control moves (
19). An advanced trajectory linearisation method adopted from [
35] is used to improve the derivation of a simple quadratic optimisation task from the nonlinear MPC-L
optimisation task.
4.1. Analytical Approximation of the MPC-L Cost-Function
In this work, it is postulated that the absolute value of the predicted error is replaced by a square of an auxiliary function
over the whole prediction horizon, i.e., for all
. Two versions of analytical approximators of the absolute value function are considered. The first one (v. 1) is described by the following equation [
36]
The second one (v. 2) is characterised by Equation [
37]
A test of both analytical approximators is necessary to find the appropriate value of parameter
c which is used in Equations (
30) and (
31).
Figure 2 depicts the comparison between the real absolute function and two approximators for different example values of the parameter
c. For the first approximator, two values are chosen for further comparisons of MPC algorithms:
and
. The value
gives quite a rough approximation; it will be interesting to analyse how this inaccuracy influences the resulting control quality. The influence of this imprecise approximation on control quality possible in MPC is discussed in
Section 5. The value
gives very good accuracy. The value
is chosen for the second approximator as it gives good accuracy.
For the first version of the approximator, the MPC-L
cost function (MPC-L
v. 1) is given by
and for the second version, the MPC-L
cost function (MPC-L
v. 2) is
4.2. Advanced Trajectory Linearisation of the MPC-L Cost-Function
As a result of using differentiable approximators of the absolute value function, both
cost-functions with the analytical approximators, as defined by Equations (
32) and (
33), are differentiable now, but still nonlinear. The source of nonlinearity is the model of the controlled process. Namely, the predictions of the future values of the controlled variable and hence the predicted control errors are nonlinear in terms of the calculated future control increments (
19). In order to solve this problem, an advanced on-line trajectory linearisation is carried out. It may be proved elsewhere [
30] that the linear approximation of the predicted trajectory of the control error over the prediction horizon embedded in the nonlinear function
can be written in general form
The future trajectory of the control error
and the trajectory of controller errors embedded in the function
are both vectors of length
N. The control errors used in the trajectory (
35) are calculated from Equation (
22) while the following equation
is utilised to determine the control errors in the trajectory (
36). Linearisation is carried out at each discrete time sampling instant
k about the following trajectory of future control scenario
which is a vector of length
. Matrix
is of size
and the vector
is of length
. The
matrix
is required for trajectory approximation (Equation (
34)). Its entries, i.e., the partial derivatives
are calculated for all
and
at each sampling instant. Differentiation depends on the type of the absolute value function approximator used. For the first approximation method (Equation (
30)), we have
and for the second one (Equation (
31)), we derive
Finally, the partial derivatives
are calculated for all
and
at each sampling instant for a given dynamical model of the process used for prediction in MPC. As far as the neural Wiener model of the PEM fuel cell system described in
Section 2 is concerned, all implementation details are given in [
24].
4.3. Formulation of the Computationally Simple MPC-L Quadratic Optimisation Task
We consider the general nonlinear MPC-L
optimisation task defined by Equation (
28), where the first part of the minimised cost function is approximated by Equation (
34). The same approximation of the trajectory
is used for both analytical approximators of the absolute value function described by Equations (
30) and (
31). Thanks to the fact that the trajectory (
34) is linear in terms of the calculated decision variable vector of MPC,
, we derive the following quadratic programming MPC-L
optimisation problem
The matrix
is of size
. Due to linearisation, the minimised cost function is quadratic in terms of the decision vector,
, and all constraints are linear with respect to the vector
. The vector constraints are defined by the following vectors of length
We name the obtained control scheme the MPC algorithm with Nonlinear Prediction and Linearisation along the Trajectory (MPC-NPLT).
The presented MPC-NPLT-L
algorithm is very universal since different kinds of dynamical models may be used for prediction. A specific form of a neural Wiener model (as shown in
Figure 1) is used in this work. Nevertheless, due to a universal formulation of the MPC-NPLT-L
algorithm, neural networks and fuzzy models of different structures may be used. The general calculation scheme, including the formulation of the quadratic optimisation task (
43), is unchanged. The actual model structure must be taken into account when the following two components are calculated:
- (a)
The predicted trajectory
for
in Equation (
37);
- (b)
The derivatives of the predicted trajectory
and
in Equation (
41) or (
42).
5. Simulations
5.1. Compared Algorithms
The effectiveness of several MPC algorithms for the PEM fuel cell system has been verified using MATLAB software. All experiments have been conducted for
,
and
; the values of MPC tuning parameters are the same as those used in previous research [
24]. Let us remind that the first version of the approximation is relatively imprecise while the second one is excellent, as depicted in
Figure 2. The MPC algorithms have been assigned to two groups.
The first group of compared control schemes consists of MPC algorithms with Nonlinear Optimisation (MPC-NO) repeated at each sampling instant. The following acronyms are used:
MPC-NO-L refers to the MPC algorithm with the classical L cost-function that measures the sum of squared predicted control errors.
MPC-NO-L refers to the MPC algorithm with the sum of absolute values of predicted control errors. It is important to stress that the real absolute function is used, without any approximation. The absolute value function is not differentiable at 0. It may create problems for a nonlinear solver used to minimise the MPC-NO-L optimisation task.
MPC-NO-L
v. 1 refers to the MPC algorithm in which the possibly non-differentiable absolute value function is replaced by its approximation defined by Equation (
30).
Finally, MPC-NO-L
v. 2 refers to the MPC algorithm in which the absolute value function is replaced by its approximation defined by Equation (
31).
The second group of compared control schemes consists of different configurations of the MPC-NPLT algorithm in which quadratic optimisation is used at each sampling instant in place of nonlinear programming necessary in the MPC-NO approach. The following acronyms are used:
As recommended in [
24], linearisation is performed around the trajectory
(Equation (
38)) defined by the most recent value of the manipulated variable computed and applied to the process at the previous sampling step, i.e.,
.
Three experiments have been conducted for the same control scenario. In the scenario, the set-point trajectory is constant, but there are five step-changes in the disturbance signal. The disturbance trajectory is presented in
Figure 3.
All simulation results of MPC algorithms which are presented in the following figures, use two analytical approximations of the absolute value function, defined by Equations (
30) and (
30), respectively, are developed for the parameter
. Additionally, a different value is considered in
Section 5.5 that compares MPC algorithms in terms of different control quality criteria.
5.2. Experiment 1
The first experiment compares the algorithms MPC-NO-L, MPC-NO-L v. 1 and MPC-NO-L v. 2. This experiment aims to check if the algorithms with approximations of the absolute value function give similar control quality to those possible when the algorithm with the classical form of the absolute value function is used. All considered MPC algorithms use nonlinear optimisation. For this purpose, the fmincon function available in MATLAB is used.
All obtained simulation results are presented in
Figure 4. The top panel presents the constant output set-point, denoted as
, and the process output signals and the middle panel depicts the trajectories of the manipulated variable. Two bottom panels depict enlarged fragments of the first panel. It is clearly shown that the recorded trajectories are nearly identical. It is evident taking into account two bottom plots, which show parts of the trajectories for two selected step changes of the disturbance. Although the difference between them is minimal, we may note that the second version of the absolute value function approximator used in the MPC-NO-L
v. 2 algorithm gives the same trajectories as those possible in the rudimentary MPC-NO-L
algorithm. The first version of the absolute value function approximator used in the MPC-NO-L
v. 1 algorithm gives slightly different trajectories. It is because the first version of the approximator gives worse results than the second one, as depicted in
Figure 2. The conclusion is that the approximators of the absolute value function give good results. Of course, nonlinear optimisation is used. This difficulty is removed in the following experiments.
5.3. Experiment 2
The following algorithms are compared in the second experiment: MPC-NO-L, MPC-NPLT-L v. 1 and MPC-NPLT-L v. 2. This experiment aims to check if the algorithms with advanced linearisation of the predicted output trajectory perform similarly to the MPC-NO-L control scheme. It is to be noted that this is the first test of the algorithms with quadratic optimisation tasks against the algorithm with the nonlinear optimisation task. Both versions of the MPC-NPLT-L algorithms use quadratic optimisation. For this purpose, the quadprog function available in MATLAB is used.
The obtained trajectories are shown in
Figure 5. The presented trajectories have minor differences between them. We can draw two conclusions. Firstly, the trajectories of the MPC-NPLT-L
v. 2 algorithm are practically the same as those generated by the MPC-NO-L
method. It means that advanced trajectory linearisation and quadratic programming may be used to eliminate nonlinear optimisation used at each sampling instant. Secondly, the trajectories of the MPC-NPLT-L
v. 1 algorithm are slightly different from two other tested algorithms. As in the first experiment, it is because the first version of the approximator gives worse results than the second one, as depicted in
Figure 2.
5.4. Experiment 3
The third experiment compares process trajectories for MPC-NPLT-L, MPC-NPLT-L v. 1 and MPC-NPLT-L v. 2 control methods. The focus is to check the performance of computationally efficient algorithms with the classical cost-function (the MPC-L cost-function) vs. two versions of approximations of the absolute value functions. All algorithms use quadratic optimisation. For this purpose, the quadprog function is used.
The trajectories of all MPC algorithms are presented in
Figure 6. We can make two observations. Firstly, in general, the algorithms with the L
cost-function give better results than the algorithm with the L
cost-function. The difference in control quality is quite significant. Secondly, as far as the algorithms with the L
norm are concerned, better results are possible when the more precise second version of the approximator is used.
5.5. Control Quality Evaluation
In this work, a multi-criteria control quality assessment is considered [
38]. As many as four control quality indicators are calculated after simulation and their values are compared. Two of these indices directly correspond with the L
and L
cost-functions. Firstly, we consider the sum of absolute values of control errors over the whole simulation horizon
and the sum of squared values of control errors
The third index is the Gauss standard deviation of the control error denoted as
. Finally, the last one is the rational entropy of the control error denoted as
H. As emphasised in [
38], control quality assessment should be performed using various indicators as only one index, e.g., the commonly used sum of squared control errors, may lead to wrong conclusions about the effectiveness of the analysed control algorithm.
Table 2 and
Table 3 compare all discussed MPC algorithms in terms of four considered control quality indicators. The first version of the absolute value function approximator uses two values of the parameter
c, i.e., 0.01 and 0.0001. Let us remind that the
gives a pretty rough approximation. The second version of the approximator uses the parameter
, which has excellent accuracy. We may make the following observations:
In general, minimisation of the L cost function gives better all control quality indicators in comparison with the L cost function, i.e., the MPC-NO-L and MPC-NPLT-L algorithms lead to bigger values of quality indices than the algorithm that the L cost function.
Advanced trajectory linearisation and quadratic optimisation in the MPC-NPLT-L algorithm lead to very similar control quality indicators as the computationally difficult nonlinear optimisation used in the MPC-NO-L scheme.
The first version of the approximator of the absolute value function with the parameter
gives worse results than the same approximator with the parameter
and the second one with
. It is straightforward since the first approximator is quite imprecise for
, as depicted in
Figure 2.
As a result of excellent accuracy, the first version of the approximator of the absolute value function with the parameter and the second version of the approximator with give very good results, practically the same.
The first version of the approximator seems to be computationally simpler. Hence, it is recommended to use the first approximator.
6. Conclusions
This paper describes a fast MPC algorithm with the L norm and its application to control the PEM fuel cell system. The classical non-differentiable absolute value function is replaced with its analytical approximation. Moreover, linearisation of the predicted process trajectory is performed at every sampling time instant to improve computational efficiency further. Two approximation methods are described and implementation details are derived. The obtained trajectories and four control quality indices have been compared to analyse the algorithm’s performance thoroughly.
The presented MPC-NPLT-L
algorithm allows to control the PEM fuel cell effectively. Firstly, it gives excellent control quality, the same as the MPC-NO-L
scheme with nonlinear optimisation. Secondly, minimisation of the L
cost function gives better results than the use of the classical L
cost function. Thirdly, the use of the analytical approximation of the absolute value function is much easier than the application of a neural approximation considered in the literature [
30]. Two versions of the analytical approximators are considered in this work. When properly tuned, both give excellent results, but the first one is computationally simpler and recommended.
It is necessary to emphasise a universal nature of the presented MPC-NPLT-L algorithm since different kinds of dynamical models may be used for prediction. Although a specific form of a neural Wiener model is used in this work, neural networks and fuzzy models of different structures may also be used. The general calculation scheme, including the formulation of the quadratic optimisation task, is unchanged. The actual model structure must be taken into account when the predicted output trajectory and the derivatives of that trajectory are calculated. The only limitation is differentiability of the model.
In the future, it is planned to study the choice of alternative cost-functions to control PEM fuel cell systems and their influence on control quality.