# Nonlinear Recurrent Neural Network Predictive Control for Energy Distribution of a Fuel Cell Powered Robot

Authors: Qihong Chen; Rong Long; Shuhai Quan; Liyan Zhang
Journal: The Scientific World Journal
 (2014)

---

## Abstract

This paper presents a neural network predictive control strategy to optimize power distribution for a fuel cell/ultracapacitor hybrid power system of a robot. We model the nonlinear power system by employing time variant auto-regressive moving average with exogenous (ARMAX), and using recurrent neural network to represent the complicated coefficients of the ARMAX model. Because the dynamic of the system is viewed as operating- state- dependent time varying local linear behavior in this frame, a linear constrained model predictive control algorithm is developed to optimize the power splitting between the fuel cell and ultracapacitor. The proposed algorithm significantly simplifies implementation of the controller and can handle multiple constraints, such as limiting substantial fluctuation of fuel cell current. Experiment and simulation results demonstrate that the control strategy can optimally split power between the fuel cell and ultracapacitor, limit the change rate of the fuel cell current, and so as to extend the lifetime of the fuel cell.


---

## Body

## 1. Introduction

As the rapid development of modern industrial technology, Ocean technology, and space technology, more and more mobile robots are demanded in these areas. Because of the advantages in operating time, weight, and dimensions, proton exchange membrane (PEM) fuel cells have been considered as alternative power sources for mobile robots.A mobile robot usually has multiple freedoms, which cause the electric load drastically to fluctuate. Sudden changes in power may significantly reduce the operating life of fuel cells in a long term [1, 2]. Furthermore, fuel cells have the characteristics of unidirectional power flow and they cannot absorb the energy from regenerative braking of a robot. As a result, fuel cells are usually arranged with auxiliary power sources to form hybrid power systems and drive mobile robots. Ultracapacitors are highly suitable for the bulk of the transient power demands since the charge/discharge current of an ultracapacitor can vary in a wide range. In this paper we choose a bank of ultracapacitors as auxiliary power source.A smart power split strategy is indispensable to enhance performance and lifetime of the hybrid power system. Jiang et al. [3] presented an adaptive control algorithm that adjusted the output current set point of the fuel cell. Ferreira et al. [4], Li et al. [5], and Kim et al. [6] developed a fuzzy controller to optimally distribute the power between the fuel cell and the battery. Rodatz et al. [7] designed an optimal control strategy to minimize the hydrogen consumption in a hybrid fuel cell system. Paladini et al. [8] proposed an optimal control strategy to power a vehicle with both fuel cell and battery to reduce fuel consumption. Lin et al. [9] studied a dynamic programming (DP) algorithm based on the fuel consumption and exhaust gas emission for a parallel electric vehicle. These strategies are effective in dealing with system efficiency but address little the lifetime of the fuel cell stack due to rapid load demand variations. Zhang et al. [10] presented a wavelet-transform algorithm to identify and allocate power demands with different frequency contents to corresponding sources to achieve an optimal power management control algorithm. This algorithm can protect fuel cell effectively but is complex and difficult to apply online. Xu et al. [11, 12] and Simmons et al. [13] proposed optimal real-time energy management strategies for a proton electrolyte membrane (PEM) fuel cell bus based on the Pontryagin’s Minimal Principle and the determined dynamic programming (DDP). Ziogou et al. [14] deployed a dynamic optimization approach based on nonlinear model of fuel cell. Li et al. [15] developed a constrained model predictive control of a solid oxide fuel cell based on genetic optimization.Undoubtedly, the fuel cell power systems are nonlinear. Therefore, the global optimization based energy management strategies depend on nonlinear models of the fuel cell power systems and are time costly. Model predictive control (MPC) has been recognized as a powerful methodology for controlling a wide class of nonlinear dynamic system [16]. In this paper we use MPC appropriately, distribute power between the fuel cell and ultracapacitor, avoid frequent fluctuation of fuel cell current, and so enhance the transient performance and extend the operating life of the hybrid system.There have been three main methods for nonlinear system modeling and predictive control [17]. The first one uses a piecewise linearization to describe the nonlinear behavior of a system. Each model is effective only in a small region, which results in that a mass of models is required [18]. The second one directly employs nonlinear models, but these involve a nonlinear online optimization problem with constraints, which is usually time-consuming and may even be unable to guarantee a feasible solution for real time control [19]. The third method is to use a local linearization approach representing a nonlinear plant, which is valid and simplifies the implement [20–24].This paper proposes an ARMAX (Autoregressive Moving Average with Exogenous input) modeling approach for fuel cell power systems. Time-variant coefficients of the ARMAX model are estimated by a recurrent neural network. The RNN-ARMAX model is an equal linear model of the fuel cell power system. Therefore, we design linear constrained model predictive control based on the RNN-ARMAX model for the nonlinear fuel cell power system. The design and implementation of the controller are significantly simplified and the method can protect fuel cell from substantial fluctuation of current by trading off transient current demand from the fuel cell to the ultracapacitor, according to constraints and weighting matrices of the output errors.The remainder of this paper is organized as follows. Section2 describes RNN-ARMAX modeling of the fuel cell power system. MPC is designed in Section 3. In Section 4, we implement and discuss simulation results. Conclusions are given in Section 5.


## 2. RNN-ARMAX Modeling

We aim at the optimization of electric power distribution between the fuel cell and ultracapacitor of a fuel cell robot.


### 2.1. System Structure and Description

The fuel cell power system studied in this paper, as shown in Figure1, is designed for a mobile robot. The electrical output of the PEM fuel cell is connected to the load through a unidirectional DC/DC converter, and an ultracapacitor bank is also connected to the load through a bidirectional DC/DC converter to form a hybrid fuel cell system. The ultracapacitor bank should supply peak power and be recharged by the fuel cell.Figure 1
Fuel cell power system of a robot.The distribution of power between the fuel cell and the ultracapacitor depends on the duty ratio of the DC/DC converters. Duty ratio of a DC/DC converter is defined as the ratio of switch on time interval,T
ON, to switching period T; that is,
(1)
d
=
T
ON
T
.There is one duty ratio,d
fc, in the unidirectional DC/DC converter for controlling output power of the fuel cell. In the bidirectional DC/DC converter, one duty ratio, d
c, is for charging the ultracapacitor, and the other, d
d, is for discharging the ultracapacitor. Power distribution is optimized by controlling the three duty ratios.


### 2.2. Identification

The hybrid system is a multiple input and multiple output nonlinear system. The control input variables are three duty ratios of the power converters. Input variables are expressed as(2)
u
(
t
)
=
(
d
fc
d
d
d
c
)
T
.The output variables contain output voltage of the fuel cell and the state of charge of the ultracapacitor and so forth. Output variables are chosen as(3)
y
(
t
)
=
(
V
fc
I
fc
I
c
SOC
V
b
I
b
)
T
,
where V
fc is voltage of the fuel cell, I
fc is current of the fuel cell, I
c is current of the ultracapacitor, SOC is state of charge of the ultracapacitor, V
b is the bus voltage: and I
b is the bus current, respectively. Power demanded by the load, P
d, is viewed as a disturbance to the system. We can describe the model as the following nonlinear function:
(4)
y
(
t
)
=
f
(
φ
(
t
)
)
+
ξ
(
t
)
,
(5)
φ
(
t
)
=
[
y
T
(
t
-
1
)
…
y
T
(
t
-
n
)
u
T
(
t
)
…
u
T
(
t
-
m
+
1
)
]
T
≜
[
φ
i
(
t
)
]
[
n
*
n
y
+
(
m
+
1
)
*
n
u
]
×
1
,
whereφ
(
t
)is the regression vector with known order n and m,n
yandn
uare dimensions of output and input,ξ
(
t
)is the system disturbance, andf
(
·
)is an unknown nonlinear function, respectively.If we design MPC based on direct use of the nonlinear model, it involves the online solution of a higher order nonlinear optimization problem with constraints, which is usually computationally expensive and may even be unable to guarantee a feasible solution for real time control.Here we use RNN-ARMAX to model the system. Performing Taylor expansion on the nonlinear functionf
(
φ
(
t
)
)around the regionφ
(
t
)
=
0as
(6)
y
(
t
)
=
f
(
0
)
+
f
′
(
0
)
φ
(
t
)
+
1
2
φ
T
(
t
)
f
′′
(
0
)
φ
(
t
)
+
…
+
ξ
(
t
)
.We introduce the notation(7)
y
0
=
f
(
0
)
,
(8)
Θ
(
φ
(
t
)
)
=
(
f
′
(
0
)
+
1
2
φ
T
(
t
)
f
′′
(
0
)
+
…
)
T
=
[
a
1
,
t
T
a
2
,
t
T
…
a
n
,
t
T
b
0
,
t
T
…
b
m
,
t
T
]
T
,
where Θ
∈
R
n
y
×
[
n
*
n
y
+
(
m
+
1
)
*
n
u
] and the coefficientsa
i
,
t
=
a
i
(
φ
(
t
)
)
,
b
i
,
t
=
b
i
(
φ
(
t
)
)are nonlinear function ofφ
(
t
).We have a regression form of the system described by (4) as follows:
(9)
y
(
t
)
=
y
0
+
Θ
(
φ
(
t
)
)
φ
(
t
)
+
ξ
(
t
)
.Here the parameter vectorΘ
(
φ
(
t
)
)is time variant. The recurrent neural network (RNN) that consists of feed-forward and feedback connections is well known to be capable of modeling and control nonlinear system. We use RNN to estimateΘ
(
φ
(
t
)
). The recurrent neural network modeling principle is shown in Figure 2.Figure 2
RNN modeling principle.The RNN is expressed as(10)
O
(
t
)
=
K
S
(
t
)
,
S
(
t
)
=
Γ
(
W
φ
(
t
)
+
W
-
S
(
t
-
1
)
)
,
Γ
(
x
)
=
[
σ
(
x
1
)
σ
(
x
1
)
…
σ
(
x
n
h
)
]
,
σ
(
x
j
)
=
1
1
+
e
-
x
j
,
where O
(
t
)
∈
R
n
y
[
n
*
n
y
+
(
m
+
1
)
*
n
u
]is output of the RNN andK
,
W
,
W
-are weights for the RNN among the output layer, the input layer, and the hidden layer. Define n
o, n
i, and n
h as the node amounts of the output layer, the input layer, and the hidden layer, respectively.K
,
W  and W
-are expressed as
(11)
K
=
[
k
1
k
2
…
k
n
h
]
n
o
×
n
h
,
S
(
t
)
=
[
s
i
(
t
)
]
n
h
×
1
,
W
=
[
w
l
k
]
n
h
×
n
i
,
W
-
=
[
w
-
l
k
]
n
h
×
n
h
,
n
o
=
n
y
[
n
*
n
y
+
(
m
+
1
)
n
u
]
,
n
i
=
n
*
n
y
+
(
m
+
1
)
n
u
.Then the output of the system is predicted by(12)
y
^
(
t
)
=
y
0
+
Ψ
(
t
)
O
(
t
)
,
whereΨ
(
t
)
∈
R
n
y
×
n
y
[
n
*
n
y
+
(
m
+
1
)
*
n
u
] and
(13)
Ψ
(
t
)
=
[
φ
T
(
t
)
0
0
0
0
φ
T
(
t
)
0
0
⋮
⋮
⋮
⋮
0
0
⋯
φ
T
(
t
)
]
≜
[
Ψ
1
(
t
)
Ψ
2
(
t
)
⋮
Ψ
n
y
(
t
)
]
.The performance criterionψ
(
t
)of the neural network is then defined by
(14)
ψ
(
t
)
=
1
2
(
y
(
t
)
-
y
0
-
Ψ
(
t
)
O
(
t
)
)
T
(
y
(
t
)
-
y
0
-
Ψ
(
t
)
O
(
t
)
)
,
where y
(
t
) is sampled output of the system. Therefore, the weights are adjusted to reduce the cost functionψ
(
t
)to a minimum value by the gradient descent method. The weight vectors are updated along with
(15)
K
(
t
+
1
)
=
K
(
t
)
-
η
∂
ψ
(
t
)
∂
K
,
W
(
t
+
1
)
=
W
(
t
)
-
η
∂
ψ
(
t
)
∂
W
,
W
-
(
t
+
1
)
=
W
-
(
t
)
-
η
∂
ψ
(
t
)
∂
W
-
,
whereηis a positive learning rate.Letq
,
r be the quotient and remainder of i
/
[
n
*
n
y
+
(
m
+
1
)
n
u
], respectively. If h
=
0, then set r
=
[
n
*
n
y
+
(
m
+
1
)
n
u
]. Else set q
=
q
+
1.  ∂
ψ
(
t
)
/
∂
K, ∂
ψ
(
t
)
/
∂
W, and  ∂
ψ
(
t
)
/
∂
W
-  are then calculated as follows:
(16)
∂
ψ
(
t
)
∂
k
i
j
=
(
y
0
,
q
+
Ψ
q
(
t
)
O
(
t
)
-
y
q
(
t
)
)
s
j
(
t
)
φ
r
(
t
)
,
∂
ψ
(
t
)
∂
w
i
j
=
(
y
0
+
Ψ
(
t
)
O
(
t
)
-
y
(
t
)
)
Ψ
(
t
)
k
i
(
t
)
H
i
(
t
)
φ
j
(
t
)
,
∂
ψ
(
t
)
∂
w
-
i
j
=
(
y
0
+
Ψ
(
t
)
O
(
t
)
-
y
(
t
)
)
Ψ
(
t
)
k
i
(
t
)
H
i
(
t
)
s
j
(
t
-
1
)
,
where
(17)
H
i
(
t
)
=
e
-
h
i
(
t
)
(
1
+
e
-
h
i
(
t
)
)
2
,
h
i
(
t
)
=
∑
l
=
1
n
i
w
i
l
(
t
)
φ
l
(
t
)
+
∑
l
=
1
n
h
w
-
i
l
(
t
)
s
l
(
t
-
1
)
.The update rules of (15) call for a proper choice of the learning rateη. For a small value ofηthe convergence is guaranteed but the speed is slow; ifηis too big, the algorithm becomes unstable. Here we develop a guideline in selecting the learning rate properly. A discrete Lyapunov function is given by
(18)
V
(
t
)
=
1
2
e
T
(
t
)
e
(
t
)
,
where
(19)
e
(
t
)
=
y
(
t
)
-
y
0
-
Ψ
(
t
)
O
(
t
)
.Thus the change of Lyapunov function due to the training process is obtained by(20)
Δ
V
(
t
)
=
V
(
t
+
1
)
-
V
(
t
)
=
1
2
[
e
T
(
t
+
1
)
e
(
t
+
1
)
-
e
T
(
t
)
e
(
t
)
]
.The error difference due to the learning is represented by(21)
e
(
t
+
1
)
=
e
(
t
)
+
Δ
e
(
t
)
=
e
(
t
)
+
[
∂
e
(
t
)
∂
W
]
T
Δ
W
,
whereΔ
Wrepresents a change in an arbitrary weight vector.From the update rule (15),
(22)
Δ
W
=
-
η
∂
ψ
(
t
)
∂
W
=
η
∂
O
T
(
t
)
∂
W
Ψ
T
(
t
)
e
(
t
)
.Then we have the following general convergence theorem.Theorem 1.
ηis the learning rate for the weights of RNN and∥
·
∥is the usual Euclidean norm inR
n. Then the convergence is guaranteed ifηis chosen as
(23)
0
<
η
<
2
max
⁡
t
⁡
∥
Ψ
(
t
)
(
∂
O
(
t
)
/
∂
W
)
∥
2
.Proof.
From equations (20)–(22),Δ
V
(
t
)can be calculated as
(24)
Δ
V
(
t
)
=
Δ
e
T
(
t
)
[
e
(
t
)
+
1
2
Δ
e
(
t
)
]
=
[
[
∂
e
(
t
)
∂
W
]
T
η
∂
O
T
(
t
)
∂
W
Ψ
T
(
t
)
e
(
t
)
]
T
×
{
e
(
t
)
+
1
2
[
∂
e
(
t
)
∂
W
]
T
η
∂
O
T
(
t
)
∂
W
Ψ
T
(
t
)
e
(
t
)
}
=
-
e
T
(
t
)
{
η
Ψ
(
t
)
∂
O
(
t
)
∂
W
∂
O
T
(
t
)
∂
W
Ψ
T
(
t
)
-
1
2
Ψ
(
t
)
∂
O
(
t
)
∂
W
∂
O
T
(
t
)
∂
W
Ψ
T
(
t
)
*
Ψ
(
t
)
∂
O
(
t
)
∂
W
∂
O
T
(
t
)
∂
W
Ψ
T
(
t
)
}
e
(
t
)
=
-
1
2
η
[
∂
O
T
(
t
)
∂
W
Ψ
T
(
t
)
e
(
t
)
]
T
×
{
2
I
-
η
∂
O
T
(
t
)
∂
W
Ψ
T
(
t
)
Ψ
(
t
)
∂
O
(
t
)
∂
W
}
×
[
∂
O
T
(
t
)
∂
W
Ψ
T
(
t
)
e
(
t
)
]
.
To guaranteeΔ
V
(
t
)
<
0,ηshould satisfy the following inequality
(25)
2
I
-
η
∂
O
T
(
t
)
∂
W
Ψ
T
(
t
)
Ψ
(
t
)
∂
O
(
t
)
∂
W
>
0
,
(26)
η
>
0
.
From inequalities (25) and (26), we obtain
(27)
0
<
η
max
⁡
t
⁡
∥
Ψ
(
t
)
∂
O
(
t
)
∂
W
∥
2
<
2
.
Namely,ηsatisfies
(28)
0
<
η
<
2
max
⁡
t
⁡
∥
Ψ
(
t
)
(
∂
O
(
t
)
/
∂
W
)
∥
2
.
This proves the theorem.We can establish a state space model from the matrix polynomials (7), (8), and (9) by defining a state vector given by
(29)
x
(
t
)
=
[
x
1
,
t
T
x
2
,
t
T
…
x
n
,
t
T
]
T
,
x
1
,
t
=
y
(
t
)
,
x
k
,
t
=
∑
i
=
1
n
+
1
-
k
a
^
i
+
k
-
1
,
t
-
1
y
(
t
-
i
)
+
∑
i
=
1
n
+
1
-
k
b
^
i
+
k
-
1
,
t
-
1
u
(
t
-
i
)
,
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
k
=
2
,
3
,
…
,
n
.A state space model can then be given by(30)
x
(
t
+
1
)
=
A
t
x
(
t
)
+
B
t
u
(
t
)
+
Ξ
(
t
+
1
)
,
y
(
t
)
=
C
x
(
t
)
,
where
(31)
i
i
i
i
i
i
i
i
A
t
=
[
a
^
1
,
t
1
0
⋯
0
a
^
2
,
t
0
1
⋯
0
⋮
⋮
⋮
⋱
⋮
a
^
n
-
1
,
t
0
0
⋯
1
a
^
n
,
t
0
0
⋯
0
]
,
B
t
=
[
b
1
,
t
b
2
,
t
⋮
b
n
,
t
]
,
Ξ
(
t
+
1
)
=
[
y
0
+
ξ
(
t
+
1
)
0
⋮
0
]
,
C
=
[
1
0
⋮
0
]
T
.Model (30) is a state space representation of MIMO RNN-ARX model (4). The parameters inA
t andB
tare estimated by the RNN, and the state  x
(
t
) at time t can be easily obtained by (29) according to the present output y
(
t
), the past input/output data, and output of the RNN.


## 2.1. System Structure and Description

The fuel cell power system studied in this paper, as shown in Figure1, is designed for a mobile robot. The electrical output of the PEM fuel cell is connected to the load through a unidirectional DC/DC converter, and an ultracapacitor bank is also connected to the load through a bidirectional DC/DC converter to form a hybrid fuel cell system. The ultracapacitor bank should supply peak power and be recharged by the fuel cell.Figure 1
Fuel cell power system of a robot.The distribution of power between the fuel cell and the ultracapacitor depends on the duty ratio of the DC/DC converters. Duty ratio of a DC/DC converter is defined as the ratio of switch on time interval,T
ON, to switching period T; that is,
(1)
d
=
T
ON
T
.There is one duty ratio,d
fc, in the unidirectional DC/DC converter for controlling output power of the fuel cell. In the bidirectional DC/DC converter, one duty ratio, d
c, is for charging the ultracapacitor, and the other, d
d, is for discharging the ultracapacitor. Power distribution is optimized by controlling the three duty ratios.


## 2.2. Identification

The hybrid system is a multiple input and multiple output nonlinear system. The control input variables are three duty ratios of the power converters. Input variables are expressed as(2)
u
(
t
)
=
(
d
fc
d
d
d
c
)
T
.The output variables contain output voltage of the fuel cell and the state of charge of the ultracapacitor and so forth. Output variables are chosen as(3)
y
(
t
)
=
(
V
fc
I
fc
I
c
SOC
V
b
I
b
)
T
,
where V
fc is voltage of the fuel cell, I
fc is current of the fuel cell, I
c is current of the ultracapacitor, SOC is state of charge of the ultracapacitor, V
b is the bus voltage: and I
b is the bus current, respectively. Power demanded by the load, P
d, is viewed as a disturbance to the system. We can describe the model as the following nonlinear function:
(4)
y
(
t
)
=
f
(
φ
(
t
)
)
+
ξ
(
t
)
,
(5)
φ
(
t
)
=
[
y
T
(
t
-
1
)
…
y
T
(
t
-
n
)
u
T
(
t
)
…
u
T
(
t
-
m
+
1
)
]
T
≜
[
φ
i
(
t
)
]
[
n
*
n
y
+
(
m
+
1
)
*
n
u
]
×
1
,
whereφ
(
t
)is the regression vector with known order n and m,n
yandn
uare dimensions of output and input,ξ
(
t
)is the system disturbance, andf
(
·
)is an unknown nonlinear function, respectively.If we design MPC based on direct use of the nonlinear model, it involves the online solution of a higher order nonlinear optimization problem with constraints, which is usually computationally expensive and may even be unable to guarantee a feasible solution for real time control.Here we use RNN-ARMAX to model the system. Performing Taylor expansion on the nonlinear functionf
(
φ
(
t
)
)around the regionφ
(
t
)
=
0as
(6)
y
(
t
)
=
f
(
0
)
+
f
′
(
0
)
φ
(
t
)
+
1
2
φ
T
(
t
)
f
′′
(
0
)
φ
(
t
)
+
…
+
ξ
(
t
)
.We introduce the notation(7)
y
0
=
f
(
0
)
,
(8)
Θ
(
φ
(
t
)
)
=
(
f
′
(
0
)
+
1
2
φ
T
(
t
)
f
′′
(
0
)
+
…
)
T
=
[
a
1
,
t
T
a
2
,
t
T
…
a
n
,
t
T
b
0
,
t
T
…
b
m
,
t
T
]
T
,
where Θ
∈
R
n
y
×
[
n
*
n
y
+
(
m
+
1
)
*
n
u
] and the coefficientsa
i
,
t
=
a
i
(
φ
(
t
)
)
,
b
i
,
t
=
b
i
(
φ
(
t
)
)are nonlinear function ofφ
(
t
).We have a regression form of the system described by (4) as follows:
(9)
y
(
t
)
=
y
0
+
Θ
(
φ
(
t
)
)
φ
(
t
)
+
ξ
(
t
)
.Here the parameter vectorΘ
(
φ
(
t
)
)is time variant. The recurrent neural network (RNN) that consists of feed-forward and feedback connections is well known to be capable of modeling and control nonlinear system. We use RNN to estimateΘ
(
φ
(
t
)
). The recurrent neural network modeling principle is shown in Figure 2.Figure 2
RNN modeling principle.The RNN is expressed as(10)
O
(
t
)
=
K
S
(
t
)
,
S
(
t
)
=
Γ
(
W
φ
(
t
)
+
W
-
S
(
t
-
1
)
)
,
Γ
(
x
)
=
[
σ
(
x
1
)
σ
(
x
1
)
…
σ
(
x
n
h
)
]
,
σ
(
x
j
)
=
1
1
+
e
-
x
j
,
where O
(
t
)
∈
R
n
y
[
n
*
n
y
+
(
m
+
1
)
*
n
u
]is output of the RNN andK
,
W
,
W
-are weights for the RNN among the output layer, the input layer, and the hidden layer. Define n
o, n
i, and n
h as the node amounts of the output layer, the input layer, and the hidden layer, respectively.K
,
W  and W
-are expressed as
(11)
K
=
[
k
1
k
2
…
k
n
h
]
n
o
×
n
h
,
S
(
t
)
=
[
s
i
(
t
)
]
n
h
×
1
,
W
=
[
w
l
k
]
n
h
×
n
i
,
W
-
=
[
w
-
l
k
]
n
h
×
n
h
,
n
o
=
n
y
[
n
*
n
y
+
(
m
+
1
)
n
u
]
,
n
i
=
n
*
n
y
+
(
m
+
1
)
n
u
.Then the output of the system is predicted by(12)
y
^
(
t
)
=
y
0
+
Ψ
(
t
)
O
(
t
)
,
whereΨ
(
t
)
∈
R
n
y
×
n
y
[
n
*
n
y
+
(
m
+
1
)
*
n
u
] and
(13)
Ψ
(
t
)
=
[
φ
T
(
t
)
0
0
0
0
φ
T
(
t
)
0
0
⋮
⋮
⋮
⋮
0
0
⋯
φ
T
(
t
)
]
≜
[
Ψ
1
(
t
)
Ψ
2
(
t
)
⋮
Ψ
n
y
(
t
)
]
.The performance criterionψ
(
t
)of the neural network is then defined by
(14)
ψ
(
t
)
=
1
2
(
y
(
t
)
-
y
0
-
Ψ
(
t
)
O
(
t
)
)
T
(
y
(
t
)
-
y
0
-
Ψ
(
t
)
O
(
t
)
)
,
where y
(
t
) is sampled output of the system. Therefore, the weights are adjusted to reduce the cost functionψ
(
t
)to a minimum value by the gradient descent method. The weight vectors are updated along with
(15)
K
(
t
+
1
)
=
K
(
t
)
-
η
∂
ψ
(
t
)
∂
K
,
W
(
t
+
1
)
=
W
(
t
)
-
η
∂
ψ
(
t
)
∂
W
,
W
-
(
t
+
1
)
=
W
-
(
t
)
-
η
∂
ψ
(
t
)
∂
W
-
,
whereηis a positive learning rate.Letq
,
r be the quotient and remainder of i
/
[
n
*
n
y
+
(
m
+
1
)
n
u
], respectively. If h
=
0, then set r
=
[
n
*
n
y
+
(
m
+
1
)
n
u
]. Else set q
=
q
+
1.  ∂
ψ
(
t
)
/
∂
K, ∂
ψ
(
t
)
/
∂
W, and  ∂
ψ
(
t
)
/
∂
W
-  are then calculated as follows:
(16)
∂
ψ
(
t
)
∂
k
i
j
=
(
y
0
,
q
+
Ψ
q
(
t
)
O
(
t
)
-
y
q
(
t
)
)
s
j
(
t
)
φ
r
(
t
)
,
∂
ψ
(
t
)
∂
w
i
j
=
(
y
0
+
Ψ
(
t
)
O
(
t
)
-
y
(
t
)
)
Ψ
(
t
)
k
i
(
t
)
H
i
(
t
)
φ
j
(
t
)
,
∂
ψ
(
t
)
∂
w
-
i
j
=
(
y
0
+
Ψ
(
t
)
O
(
t
)
-
y
(
t
)
)
Ψ
(
t
)
k
i
(
t
)
H
i
(
t
)
s
j
(
t
-
1
)
,
where
(17)
H
i
(
t
)
=
e
-
h
i
(
t
)
(
1
+
e
-
h
i
(
t
)
)
2
,
h
i
(
t
)
=
∑
l
=
1
n
i
w
i
l
(
t
)
φ
l
(
t
)
+
∑
l
=
1
n
h
w
-
i
l
(
t
)
s
l
(
t
-
1
)
.The update rules of (15) call for a proper choice of the learning rateη. For a small value ofηthe convergence is guaranteed but the speed is slow; ifηis too big, the algorithm becomes unstable. Here we develop a guideline in selecting the learning rate properly. A discrete Lyapunov function is given by
(18)
V
(
t
)
=
1
2
e
T
(
t
)
e
(
t
)
,
where
(19)
e
(
t
)
=
y
(
t
)
-
y
0
-
Ψ
(
t
)
O
(
t
)
.Thus the change of Lyapunov function due to the training process is obtained by(20)
Δ
V
(
t
)
=
V
(
t
+
1
)
-
V
(
t
)
=
1
2
[
e
T
(
t
+
1
)
e
(
t
+
1
)
-
e
T
(
t
)
e
(
t
)
]
.The error difference due to the learning is represented by(21)
e
(
t
+
1
)
=
e
(
t
)
+
Δ
e
(
t
)
=
e
(
t
)
+
[
∂
e
(
t
)
∂
W
]
T
Δ
W
,
whereΔ
Wrepresents a change in an arbitrary weight vector.From the update rule (15),
(22)
Δ
W
=
-
η
∂
ψ
(
t
)
∂
W
=
η
∂
O
T
(
t
)
∂
W
Ψ
T
(
t
)
e
(
t
)
.Then we have the following general convergence theorem.Theorem 1.
ηis the learning rate for the weights of RNN and∥
·
∥is the usual Euclidean norm inR
n. Then the convergence is guaranteed ifηis chosen as
(23)
0
<
η
<
2
max
⁡
t
⁡
∥
Ψ
(
t
)
(
∂
O
(
t
)
/
∂
W
)
∥
2
.Proof.
From equations (20)–(22),Δ
V
(
t
)can be calculated as
(24)
Δ
V
(
t
)
=
Δ
e
T
(
t
)
[
e
(
t
)
+
1
2
Δ
e
(
t
)
]
=
[
[
∂
e
(
t
)
∂
W
]
T
η
∂
O
T
(
t
)
∂
W
Ψ
T
(
t
)
e
(
t
)
]
T
×
{
e
(
t
)
+
1
2
[
∂
e
(
t
)
∂
W
]
T
η
∂
O
T
(
t
)
∂
W
Ψ
T
(
t
)
e
(
t
)
}
=
-
e
T
(
t
)
{
η
Ψ
(
t
)
∂
O
(
t
)
∂
W
∂
O
T
(
t
)
∂
W
Ψ
T
(
t
)
-
1
2
Ψ
(
t
)
∂
O
(
t
)
∂
W
∂
O
T
(
t
)
∂
W
Ψ
T
(
t
)
*
Ψ
(
t
)
∂
O
(
t
)
∂
W
∂
O
T
(
t
)
∂
W
Ψ
T
(
t
)
}
e
(
t
)
=
-
1
2
η
[
∂
O
T
(
t
)
∂
W
Ψ
T
(
t
)
e
(
t
)
]
T
×
{
2
I
-
η
∂
O
T
(
t
)
∂
W
Ψ
T
(
t
)
Ψ
(
t
)
∂
O
(
t
)
∂
W
}
×
[
∂
O
T
(
t
)
∂
W
Ψ
T
(
t
)
e
(
t
)
]
.
To guaranteeΔ
V
(
t
)
<
0,ηshould satisfy the following inequality
(25)
2
I
-
η
∂
O
T
(
t
)
∂
W
Ψ
T
(
t
)
Ψ
(
t
)
∂
O
(
t
)
∂
W
>
0
,
(26)
η
>
0
.
From inequalities (25) and (26), we obtain
(27)
0
<
η
max
⁡
t
⁡
∥
Ψ
(
t
)
∂
O
(
t
)
∂
W
∥
2
<
2
.
Namely,ηsatisfies
(28)
0
<
η
<
2
max
⁡
t
⁡
∥
Ψ
(
t
)
(
∂
O
(
t
)
/
∂
W
)
∥
2
.
This proves the theorem.We can establish a state space model from the matrix polynomials (7), (8), and (9) by defining a state vector given by
(29)
x
(
t
)
=
[
x
1
,
t
T
x
2
,
t
T
…
x
n
,
t
T
]
T
,
x
1
,
t
=
y
(
t
)
,
x
k
,
t
=
∑
i
=
1
n
+
1
-
k
a
^
i
+
k
-
1
,
t
-
1
y
(
t
-
i
)
+
∑
i
=
1
n
+
1
-
k
b
^
i
+
k
-
1
,
t
-
1
u
(
t
-
i
)
,
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
k
=
2
,
3
,
…
,
n
.A state space model can then be given by(30)
x
(
t
+
1
)
=
A
t
x
(
t
)
+
B
t
u
(
t
)
+
Ξ
(
t
+
1
)
,
y
(
t
)
=
C
x
(
t
)
,
where
(31)
i
i
i
i
i
i
i
i
A
t
=
[
a
^
1
,
t
1
0
⋯
0
a
^
2
,
t
0
1
⋯
0
⋮
⋮
⋮
⋱
⋮
a
^
n
-
1
,
t
0
0
⋯
1
a
^
n
,
t
0
0
⋯
0
]
,
B
t
=
[
b
1
,
t
b
2
,
t
⋮
b
n
,
t
]
,
Ξ
(
t
+
1
)
=
[
y
0
+
ξ
(
t
+
1
)
0
⋮
0
]
,
C
=
[
1
0
⋮
0
]
T
.Model (30) is a state space representation of MIMO RNN-ARX model (4). The parameters inA
t andB
tare estimated by the RNN, and the state  x
(
t
) at time t can be easily obtained by (29) according to the present output y
(
t
), the past input/output data, and output of the RNN.


## 3. Controller Design

A predictive controller will be designed to predict the output trajectory of the fuel cell power system and compute a series of control actions, subject to constraints, that will minimize the difference between the predicted trajectory and desired trajectory. A prominent advantage of this controller over other control schemes is its ability to deal with constraints in a systematic and straightforward manner.To design predictive controller for the system, an objective function is defined as [18]
(32)
min
⁡
u
(
t
)
,
⋯
,
u
(
t
+
N
-
1
)
⁡
J
=
∑
k
=
1
N
u
(
(
y
^
(
t
+
k
)
-
y
r
(
t
+
k
)
)
T
i
i
i
i
i
×
Q
(
y
^
(
t
+
k
)
-
y
r
(
t
+
k
)
)
iiii
i
+
u
T
(
t
+
k
)
R
u
T
(
t
+
k
)
)
,
where N
u is predictive horizon, y
^
(
t
+
k
) is the estimated output of the system at instant t
+
k through models based on information available at instantt. y
r
(
t
+
k
) is the desired output at instant t
+
k, and Q
,
R are weighting matrices on output errors and control, respectively. We choose the control horizon to be equal to the prediction horizon and defineQ
=
diag
⁡
(
Q
V
fc
Q
I
fc
Q
I
c
Q
SOC
Q
V
b
)andR
=
diag
⁡
(
R
d
fc
R
d
d
R
d
c
), where,Q
V
fc
,
Q
I
fc
,
Q
I
c
,
Q
SOC, and Q
V
b are penalties on errors in  V
fc
,
I
fc
,
I
c
,
SOC andV
b, respectively.R
d
fc
,
R
d
d, andR
d
care penalties ond
fc
,
d
d and d
c, respectively.Substituting state equations (30) into (32), the equation is abbreviated as
(33)
min
⁡
U
⁡
J
=
1
2
∥
U
∥
Ω
2
+
(
y
T
(
t
)
L
+
G
)
U
,
where y
(
t
) is system output at instant t, and Ω
,
L
,
G are constant matrices calculated through the system model and matrices Q
,
R.Consider the following:(34)
U
=
(
u
T
(
t
)
,
u
T
(
t
+
1
)
,
…
,
u
T
(
t
+
N
u
-
1
)
)
.In the hybrid system, there are several limits to deal with. Rapid variation on current will reduce lifetime of fuel cell, so it is required to constrain the fluctuation of fuel cell current; that is,(35)
-
Δ
I
max
⁡
≤
I
fc
(
t
+
1
)
-
I
fc
(
t
)
≤
Δ
I
max
⁡
,
where Δ
I
max
⁡ is the acceptable maximum value. Moreover, the state of charge of the ultracapacitor, the current of the ultracapacitor, and the voltage of the fuel cell should be limited to some expected range:
(36)
SO
C
min
⁡
≤
SOC
≤
SO
C
max
⁡
,
(37)
-
I
c
,
max
⁡
≤
I
c
≤
I
c
,
max
⁡
,
(38)
V
fc
,
min
⁡
≤
V
fc
≤
V
fc
,
max
⁡
,
where SOC
min
⁡ and V
fc
,
min
⁡ are the lower limitations, SO
C
max
⁡, I
c
,
max
⁡, and V
fc
,
max
⁡ are the upper limits, respectively. These limitations are determined by the characteristics of the ultracapacitor and fuel cell.A prominent advantage of MPC is its ability to deal with constraints. Deduced from equations (30), (32) and inequalities (35)–(38), the control optimization is transformed to the following constrained quadratic programming problem:
(39)
min
⁡
U
J
=
1
2
∥
U
∥
Ω
2
+
(
y
T
(
t
)
L
+
G
)
U
,
s
.
t
.
U
min
⁡
≤
E
U
≤
U
max
⁡
,
where  U
min
⁡
,
U
max
⁡
∈
R
m, and E
∈
R
m
×
N
u are constant matrices obtained from (30) and inequalities (35)–(38). We can solve this optimal problem using the neural network method investigated in [25].


## 4. Experiment and Simulation

The hybrid fuel cell system, as shown in Figure1, is designed to power a robot. The rated power is 500 W. The DC bus voltage is controlled around 24 V. The PEM fuel cells have 40 cells and an active area of 22 cm2. The ultracapacitor is 200 F and the rated voltage is 24 V. The value of capacitance can be realized by a bank of 8 ultracapacitors, each with capacitance of 1600 F and a rated voltage of 3 V, connected in series. The upper and lower limits of SOC are 1 and 0.45, respectively. The maximum stored energy is 16 W h, although only 12.76 W h is available between the maximum and minimum of SOC. This 12.76 W h corresponds to an average power at 500 W for 92 seconds and that is sufficient to buffer the fuel cell from acceleration transients.


### 4.1. Modeling Experiment and Simulation

When real input and output data of the PEM fuel cell was sampled, the operating parameters are shown in Table1.Table 1
Parameters used in the experiment and simulation.
Sym.
Meaning
Value
T
st
Temperature of fuel cell
343 K
T
atm
Atmospheric temperature
295 K
P
H
2
Partial pressure of hydrogen
1.5 atm
n
Number of cells in each stack
40
A
Active area of fuel cell
22 cm2
C
Capacitance of ultracapacitor
200 F
V
c
,
max
⁡
Rated voltage of ultracapacitor
24 VThe collected data are equally divided into two groups. The first group is used for modeling and the second group is used for validating. The simulated and measured V-I characteristics curves of the fuel cell are shown in Figure3. Current of the ultracapacitor changes as Figure 4, and the simulated and measured voltage curves are shown in Figure 5. It is shown that the RNN-ARMAX model closely matches the practical fuel cell power system.Figure 3
The simulated and measured V-I characteristics curves of the fuel cell.Figure 4
Current of the ultracapacitor.Figure 5
The simulated and measured voltage of the ultracapacitor.


### 4.2. Control Simulation

Control performances of constrained and unconstrained MPCs are studied and compared to validate the proposed constrained MPC. The constraints of the constrained MPC are listed in Table2.Table 2
Constraints for the constrained MPC.
Sym.
Meaning
Lower limit
Upper limit
Δ
I
max
⁡
Rate of change of fuel cell current
−0.4 A/s
0.4 A/s
SOC
State of charge of the ultracapacitor
0.45
1
I
c
Current of the ultracapacitor
−30 A
30 A
V
st
Voltage of the fuel cell
27.5 V
40 VA typical load cycle that is used in simulation and the power profile, as shown in Figure6, is considered as the power demand.Figure 6
Power profile.The simulation results for both the unconstrained and the constrained MPC are shown in Figure7. It is shown that, there exist significant perturbations in current of fuel cell for unconstrained MPC. This phenomenon may cause oxygen starvation because the dynamic response of oxygen supply is slower, while in the case of the constrained MPC, current and voltage are much smoother.Simulation results of constrained and unconstrained MPC: (a) current of fuel cell; (b) voltage of fuel cell; (c)SOC of ultracapacitor.
(a)
(b)
(c)In the case of constrained MPC, the oscillation ofSOC of the ultracapacitor is much larger than that of the unconstrained MPC. The reason is that constrained MPC draws much more energy from the ultracapacitor to supply the peak load and so limits perturbations of the current of the fuel cell.Constraint results are shown in Figure8. It’s exciting that the maximum rate of change of the fuel cell is 0.4 A/s, the minimum voltage of the fuel cell is 27.5 V, the charge and discharge current of the ultracapacitor are no more than 30 A, and the SOC of the ultracapacitor is between 0.45 and 1. It is shown that these variables change in the desired and constrained ranges. These phenomena demonstrate that the constraints on the fuel cell power system are valid.Curves for validating of constraints: (a) change rate of fuel cell current; (b) voltage of fuel cell; (c) current of ultracapacitor; (d)SOC of ultracapacitor.
(a)
(b)
(c)
(d)The power split under the constrained MPC is shown in Figure9. We set the minimum voltage of the fuel cell as 27.5 V and the corresponding maximum power of fuel cell as 500 W. It is noticed that the fuel cell power changes in low speed and is no more than 500 W. The high frequency power demands are squeezed from the ultracapacitor. Furthermore, SOC,  I
c and other constrained variables satisfy their constraints. Consequently, the output power of the fuel cell is well controlled and it is helpful to extend the operating life of the fuel cell.Figure 9
Power distribution of the hybrid system.


## 4.1. Modeling Experiment and Simulation

When real input and output data of the PEM fuel cell was sampled, the operating parameters are shown in Table1.Table 1
Parameters used in the experiment and simulation.
Sym.
Meaning
Value
T
st
Temperature of fuel cell
343 K
T
atm
Atmospheric temperature
295 K
P
H
2
Partial pressure of hydrogen
1.5 atm
n
Number of cells in each stack
40
A
Active area of fuel cell
22 cm2
C
Capacitance of ultracapacitor
200 F
V
c
,
max
⁡
Rated voltage of ultracapacitor
24 VThe collected data are equally divided into two groups. The first group is used for modeling and the second group is used for validating. The simulated and measured V-I characteristics curves of the fuel cell are shown in Figure3. Current of the ultracapacitor changes as Figure 4, and the simulated and measured voltage curves are shown in Figure 5. It is shown that the RNN-ARMAX model closely matches the practical fuel cell power system.Figure 3
The simulated and measured V-I characteristics curves of the fuel cell.Figure 4
Current of the ultracapacitor.Figure 5
The simulated and measured voltage of the ultracapacitor.


## 4.2. Control Simulation

Control performances of constrained and unconstrained MPCs are studied and compared to validate the proposed constrained MPC. The constraints of the constrained MPC are listed in Table2.Table 2
Constraints for the constrained MPC.
Sym.
Meaning
Lower limit
Upper limit
Δ
I
max
⁡
Rate of change of fuel cell current
−0.4 A/s
0.4 A/s
SOC
State of charge of the ultracapacitor
0.45
1
I
c
Current of the ultracapacitor
−30 A
30 A
V
st
Voltage of the fuel cell
27.5 V
40 VA typical load cycle that is used in simulation and the power profile, as shown in Figure6, is considered as the power demand.Figure 6
Power profile.The simulation results for both the unconstrained and the constrained MPC are shown in Figure7. It is shown that, there exist significant perturbations in current of fuel cell for unconstrained MPC. This phenomenon may cause oxygen starvation because the dynamic response of oxygen supply is slower, while in the case of the constrained MPC, current and voltage are much smoother.Simulation results of constrained and unconstrained MPC: (a) current of fuel cell; (b) voltage of fuel cell; (c)SOC of ultracapacitor.
(a)
(b)
(c)In the case of constrained MPC, the oscillation ofSOC of the ultracapacitor is much larger than that of the unconstrained MPC. The reason is that constrained MPC draws much more energy from the ultracapacitor to supply the peak load and so limits perturbations of the current of the fuel cell.Constraint results are shown in Figure8. It’s exciting that the maximum rate of change of the fuel cell is 0.4 A/s, the minimum voltage of the fuel cell is 27.5 V, the charge and discharge current of the ultracapacitor are no more than 30 A, and the SOC of the ultracapacitor is between 0.45 and 1. It is shown that these variables change in the desired and constrained ranges. These phenomena demonstrate that the constraints on the fuel cell power system are valid.Curves for validating of constraints: (a) change rate of fuel cell current; (b) voltage of fuel cell; (c) current of ultracapacitor; (d)SOC of ultracapacitor.
(a)
(b)
(c)
(d)The power split under the constrained MPC is shown in Figure9. We set the minimum voltage of the fuel cell as 27.5 V and the corresponding maximum power of fuel cell as 500 W. It is noticed that the fuel cell power changes in low speed and is no more than 500 W. The high frequency power demands are squeezed from the ultracapacitor. Furthermore, SOC,  I
c and other constrained variables satisfy their constraints. Consequently, the output power of the fuel cell is well controlled and it is helpful to extend the operating life of the fuel cell.Figure 9
Power distribution of the hybrid system.


## 5. Conclusions

RNN-ARMAX model was established and linear constrained MPC was developed and verified for a fuel cell power system. The proposed approach, different from other approaches, models the nonlinear fuel cell power system as linear time varying system. Accordingly, linear constrained MPC can be used to globally optimize power distribution and deal with limitations. The design and implementation of the controller are significantly simplified and the method can protect fuel cell from substantial fluctuation of current by trading off transient current demand from the fuel cell to the ultracapacitor.



---

*Source: 509729-2014-02-20.xml*