# Modeling, Parameters Identification, and Control of High Pressure Fuel Cell Back-Pressure Valve

Authors: Fengxiang Chen; Ling Liu; Shiguang Liu; Tong Zhang
Journal: Mathematical Problems in Engineering
 (2014)

---

## Abstract

The reactant pressure is crucial to the efficiency and lifespan of a high pressure PEMFC engine. This paper analyses a regulated back-pressure valve (BPV) for the cathode outlet flow in a high pressure PEMFC engine, which can achieve precisely pressure control. The modeling, parameters identification, and nonlinear controller design of a BPV system are considered. The identified parameters are used in designing active disturbance rejection controller (ADRC). Simulations and extensive experiments are conducted with the xPC Target and show that the proposed controller can not only achieve good dynamic and static performance but also have strong robustness against parameters’ disturbance and external disturbance.


---

## Body

## 1. Introduction

High pressure proton exchange membrane fuel cell (PEMFC) engines offer clean and efficient energy production and are advantaged of high stack power density, quick system response, and easy humidification. One of the dominant design parameters for high pressure PEMFC system is the operating pressure [1]. In the case of stored compressed hydrogen, the pressure of the cathode air supply system becomes an important optimized parameter [2]. Compared with the low-pressure system, high operating pressure improves the performance of the fuel cell stack but consumes more power due to the flow device, such as the compressor (nearly consuming up to 20 percent of the net power [3]). In order to maintain a high efficiency during operation, it is important to control air pressure effectively. In this study, an electronic throttle (ET) body is served as a back-pressure valve (BPV); the cathode air supply system can operate in a proper pressure mainly by auto-tuning the angle of the BPV for cathode outlet flow at a given power. Therefore, the fast and accurate control of the angle of the BPV is a key issue to regulate the back pressure for a high pressure PEMFC engine.The analysis of BPV model-based control system may begin with the mathematical model of a BPV system. The BPV is different from the traditional ET in that it has to consider the influence of the flow resistance. In addition, the BPV has many unknown parameters and there is no systematic method to estimate its parameters so far. Song has estimated some parameters [4]. The simulations show that the model and real plant have unacceptable bias due to excess approximation. Grep and Lee have designed some simple algorithms [5] and modified them to identify the ET parameters in [6]. But there is no mention of the parameter verification. Based on the previous studies, several algorithms are developed to estimate the unknown parameters in this paper. Parameter identification uses simple equipment, such as current, voltage sensor, and obtained parameters are optimized through the optimization toolbox. The model and acquired parameters can be used for the model-based controller design. What is more, they can also be used to find out the optimal opening of the BPV at a given power by merging into a PEMFC model. At present, set-points are usually obtained by repeating experiments on PEMFC engine.In the aspect of control, various control design techniques have been studied in the previous publications, such as proportional-integral-derivative (PID) control [7] and sliding mode control [8]. In [9], the application of receptive field weighted regression (RFWR) method is demonstrated for the composite control of ET. In [10], an ET controller based on the idea of back-stepping method is proposed, which meets the general requirements of the state stability. Then, in [11], the back-stepping controller combined with state observer is designed for the real-time system. In [12], a variable structure concept is utilized based on the sliding mode observer technique and feedback back-stepping technique. However, the above methods have strong model dependence. Parameters-varying or external disturbances can usually deteriorate the control performance. Most importantly, the system parameters used in the simulations are not verified with actual ET model.In this paper, a detailed nonlinear dynamic model of BPV is provided and several algorithms are designed to estimate the unknown parameters. The methods used in the identification process consist of extensive measurements, experiments, and verifications. The back-pressure valve used in this paper is a commercial available ET system made by Bosch, in type of DV-E5. The rest of the paper is organized as follows: in Section2 a detailed nonlinear dynamic model of BPV is provided. Next in Section 3, the identification algorithms are explained. Finally, the active disturbance rejection controller (ADRC) is designed based on the mathematical model and the identified parameters in Section 4. After that, the model controlling effect and real-time controlling effect are presented.


## 2.  Modeling of the Back-Pressure Valve

The ET is featured with high nonlinearity. The model of ET has been approached, which includes the nonlinear return spring and friction model. Yet, the BPV model is different from the traditional ET in that it has to consider the influence of the flow resistance.The BPV is a variable flow resistance device which provides needed pressure based on the position of the valve plate. It consists of a brushed DC motor with permanent magnets linked to the valve shaft by means of a set of spur gearbox, a nonlinear return spring system, a valve plate, and two redundant potentiometers (Figure1) [13].Figure 1
The structure of back-pressure valve.


### 2.1. Mechanism Modeling

The motor voltageE
a is defined as the control input and the BPV angle θ as the output, where θ is the angle between valve plate and vertical direction of inlet flow shown in Figure 2. The angle of the valve is between 10.5° and 101.5° due to the physical limit, and the effective opening range is 0–91°.Figure 2
The angle of the BPV.The BPV is driven by a brushed DC motor. The motor model can be described by the following equations [14]:
(1)
E
a
=
R
a
i
a
+
L
a
d
i
a
d
t
+
V
b
,
(2)
V
b
=
k
b
θ
˙
m
,
(3)
T
m
=
k
t
i
a
,
where V
b is the back electromotive force, and the value of L
a is so tiny that it can be ignored. The details of other variables are defined in the Nomenclature part. The mechanical equations of the motion are described by
(4)
J
m
θ
¨
m
=
-
k
m
θ
˙
m
-
T
L
+
T
m
,
J
g
θ
¨
=
-
T
sp
-
T
tf
+
n
T
L
,
where n
=
θ
m
/
θ, and T
sp and T
tf are described in the following:
(5)
T
sp
=
k
op
(
θ
-
θ
0
)
+
k
pre
θ
0
<
θ
<
θ
op
k
cl
⁡
(
θ
-
θ
0
)
-
k
pre
θ
cl
⁡
<
θ
<
θ
0
0
θ
=
θ
0
,
T
tf
=
k
f
θ
˙
+
k
tf
sgn
⁡
θ
˙
.When considering previous equations (1)–(5), the equation of motion of the electromechanical system can be written in form:
(6)
J
θ
¨
+
n
2
k
m
θ
˙
+
T
sp
+
T
tf
=
n
k
t
i
a
,
where J is the inertial moment in normalized units and J
=
J
g
+
n
2
J
m.If we definex
1
=
θ, x
2
=
θ
˙ as state variables, u
∈
[
-
E
a
,
E
a
] as the controlled input, then the system can be described by
(7)
x
˙
1
x
˙
2
=
0
1
-
k
sp
J
-
a
x
1
x
2
+
0
n
k
t
J
R
a
u
+
0
b
(
x
2
)
+
0
1
d
,
where b
(
x
2
)
=
(
-
k
tf
sgn
⁡
(
x
2
)
+
k
sp
θ
0
-
k
pre
sgn
⁡
(
x
1
-
θ
0
)
)
/
J, a
=
(
n
2
k
m
+
k
f
+
n
2
k
b
k
t
/
R
a
)
/
J, and d stands for the summation of unmodeled dynamics and unknown disturbances, for example, the load torque of inlet flow.


### 2.2.  Fluid Modeling

This model is based on the PEMFC model proposed by Pukrushpan [15]. The mass conservation principle is used to develop the outlet manifold model. For any manifold,
(8)
d
m
d
t
=
W
in
-
W
out
,
where W
in and W
out are mass flow rates in and out of the manifold, respectively.The temperature of the stack has been well regulated to a constant by the thermal management of the PEMFC. Therefore, the change of air temperature in the outlet manifold is negligible, and thus the pressure can be determined by(9)
d
p
om
d
t
=
R
g
T
om
V
om
d
m
d
t
.The outlet mass flow is governed by a BPV. The BPV is a variable flow resistance device which provides needed pressure based on the position of the valve plate. The flow is a function of manifold pressurep
om, the valve angle θ, the diameter of the BPV D, and the downstream pressure of the manifold p
atm. Thus, the air flow through the BPV satisfies [16]
(10)
Q
N
=
28.5
A
K
V
0
D
2
p
2
Δ
p
S
,
(11)
K
V
0
=
1
-
cos
⁡
θ
e
cos
⁡
2
θ
,
(12)
Δ
p
=
p
1
-
p
2
p
2
p
2
>
0.5
p
2
p
2
p
2
≤
0.5
,
where K
V
0 is the flow coefficient, which is related to the opening of the BPV, and 0
≤
K
V
0
≤
1. D is the valve diameter. A is the structure constant, which can be set to 6.2 for our BPV [16]. Δ
p
(
kgf
/
c
m
2
) is the pressure drop after the valve (which is less than half of the inlet pressure), p
1
=
p
om, p
2
=
p
atm. When the pressure drop exceeds the limit, the flow dynamic will be choked, which may lead to the instability of the flow rate [17]. Choked phenomenon is shown in Figure 3, where Δ
p
cr is the critical pressure. When pressure difference exceeds the critical pressure, for example, Δ
p
=
Δ
p
a, the actual flow is Q
max
⁡, where Q
max
⁡
′
≫
Q
max
⁡.Figure 3
Choked phenomenon.


## 2.1. Mechanism Modeling

The motor voltageE
a is defined as the control input and the BPV angle θ as the output, where θ is the angle between valve plate and vertical direction of inlet flow shown in Figure 2. The angle of the valve is between 10.5° and 101.5° due to the physical limit, and the effective opening range is 0–91°.Figure 2
The angle of the BPV.The BPV is driven by a brushed DC motor. The motor model can be described by the following equations [14]:
(1)
E
a
=
R
a
i
a
+
L
a
d
i
a
d
t
+
V
b
,
(2)
V
b
=
k
b
θ
˙
m
,
(3)
T
m
=
k
t
i
a
,
where V
b is the back electromotive force, and the value of L
a is so tiny that it can be ignored. The details of other variables are defined in the Nomenclature part. The mechanical equations of the motion are described by
(4)
J
m
θ
¨
m
=
-
k
m
θ
˙
m
-
T
L
+
T
m
,
J
g
θ
¨
=
-
T
sp
-
T
tf
+
n
T
L
,
where n
=
θ
m
/
θ, and T
sp and T
tf are described in the following:
(5)
T
sp
=
k
op
(
θ
-
θ
0
)
+
k
pre
θ
0
<
θ
<
θ
op
k
cl
⁡
(
θ
-
θ
0
)
-
k
pre
θ
cl
⁡
<
θ
<
θ
0
0
θ
=
θ
0
,
T
tf
=
k
f
θ
˙
+
k
tf
sgn
⁡
θ
˙
.When considering previous equations (1)–(5), the equation of motion of the electromechanical system can be written in form:
(6)
J
θ
¨
+
n
2
k
m
θ
˙
+
T
sp
+
T
tf
=
n
k
t
i
a
,
where J is the inertial moment in normalized units and J
=
J
g
+
n
2
J
m.If we definex
1
=
θ, x
2
=
θ
˙ as state variables, u
∈
[
-
E
a
,
E
a
] as the controlled input, then the system can be described by
(7)
x
˙
1
x
˙
2
=
0
1
-
k
sp
J
-
a
x
1
x
2
+
0
n
k
t
J
R
a
u
+
0
b
(
x
2
)
+
0
1
d
,
where b
(
x
2
)
=
(
-
k
tf
sgn
⁡
(
x
2
)
+
k
sp
θ
0
-
k
pre
sgn
⁡
(
x
1
-
θ
0
)
)
/
J, a
=
(
n
2
k
m
+
k
f
+
n
2
k
b
k
t
/
R
a
)
/
J, and d stands for the summation of unmodeled dynamics and unknown disturbances, for example, the load torque of inlet flow.


## 2.2.  Fluid Modeling

This model is based on the PEMFC model proposed by Pukrushpan [15]. The mass conservation principle is used to develop the outlet manifold model. For any manifold,
(8)
d
m
d
t
=
W
in
-
W
out
,
where W
in and W
out are mass flow rates in and out of the manifold, respectively.The temperature of the stack has been well regulated to a constant by the thermal management of the PEMFC. Therefore, the change of air temperature in the outlet manifold is negligible, and thus the pressure can be determined by(9)
d
p
om
d
t
=
R
g
T
om
V
om
d
m
d
t
.The outlet mass flow is governed by a BPV. The BPV is a variable flow resistance device which provides needed pressure based on the position of the valve plate. The flow is a function of manifold pressurep
om, the valve angle θ, the diameter of the BPV D, and the downstream pressure of the manifold p
atm. Thus, the air flow through the BPV satisfies [16]
(10)
Q
N
=
28.5
A
K
V
0
D
2
p
2
Δ
p
S
,
(11)
K
V
0
=
1
-
cos
⁡
θ
e
cos
⁡
2
θ
,
(12)
Δ
p
=
p
1
-
p
2
p
2
p
2
>
0.5
p
2
p
2
p
2
≤
0.5
,
where K
V
0 is the flow coefficient, which is related to the opening of the BPV, and 0
≤
K
V
0
≤
1. D is the valve diameter. A is the structure constant, which can be set to 6.2 for our BPV [16]. Δ
p
(
kgf
/
c
m
2
) is the pressure drop after the valve (which is less than half of the inlet pressure), p
1
=
p
om, p
2
=
p
atm. When the pressure drop exceeds the limit, the flow dynamic will be choked, which may lead to the instability of the flow rate [17]. Choked phenomenon is shown in Figure 3, where Δ
p
cr is the critical pressure. When pressure difference exceeds the critical pressure, for example, Δ
p
=
Δ
p
a, the actual flow is Q
max
⁡, where Q
max
⁡
′
≫
Q
max
⁡.Figure 3
Choked phenomenon.


## 3. Parameter Identification

The valve parameters are usually unknown or cannot be accurately measured. Advanced nonlinear control algorithm requires modeling of the plant and good estimate of its parameters. This section provides some techniques to identify the parameters of the aforementioned model [18].


### 3.1.  xPC Target

The real-time simulation of rapid control prototyping (RCP) is that the real plant is controlled by simulating model. For our RCP, we have used xPC Target system (Figure4). Various tests are conducted to estimate the parameters on the platform.Figure 4
xPC Target platform.


### 3.2. Characteristic Curve

The plant has two potentiometers supplied by 5-volt DC. The calibration curves of Throttle Position Sensor (TPS) 1 and TPS2 are shown in Figure5. The angle of the valve is between 10.5° and 101.5° due to the physical limit. When there is no voltage applied to the driven motor, the output voltage of TPS1 is 0.83 V, so the default position θ
0 is 17.5°.Figure 5
Characteristic curve of angle sensors.


### 3.3. Identification ofR
a

There is nearly no rotation when a slowly changing voltage is applied to the motor at the almost closed position [19]. Since the motor current changes slowly, the voltage across the inductance can be neglected. Besides, the inner resistance of the storage can be ignored. Therefore, the armature circuit equations (1) and (2) can be described by
(13)
E
a
=
R
a
i
a
.Figure6 displays the captured curves of the motor voltage against the motor current. Due to the good symmetry of the H-bridge, the direction of the current and voltage is neglected. The inverse of the slope represents the armature resistance which reads R
a
=
2.0
Ω.Figure 6
Motor current versus voltage.


### 3.4. Identification ofk
b

The circuit equations (1) and (2) can be rewritten as
(14)
V
b
=
n
k
b
θ
˙
=
E
a
-
R
a
i
a
,
where n
=
20.68 is the total gear ratio of the two set reduction gear. The back-EMF voltage is generated when the valve shaft is rotating. Figure 7 shows the captured back-EMF voltage and the plate’s angular velocity.Figure 7
Velocity versus back-EMF voltage.We get the back electromotive force constantk
b
=
0.0184 V/(rad/s) at the time when the plate’s angular velocity is steady. Usually, k
b and k
t are almost the same for one motor.


### 3.5.  Static Load Test

Figure8 shows the typical system response when inputting sinusoidal signal. If the frequency is low, the angular velocity and angular acceleration can be ignored. Therefore, the dynamic equation of the system can be described by
(15)
k
sp
θ
-
θ
0
+
k
pre
sgn
⁡
θ
-
θ
0
+
k
tf
sgn
⁡
θ
˙
=
n
k
t
i
a
=
T
.Figure 8
Static loads test.The response plotted in load torque versus angle characteristic is shown in Figure9. This experiment identifies the parameters: k
sp, k
pre, and k
tf.(1)
The inverse slope of the inserted straight lines represented the spring stiffness which readsk
op
=
0.4387
N
·
m
/
rad, k
cl
⁡
=
7.3693
N
·
m
/
rad.
(2)
The friction torque is half of the deviation between two inserted straight lines during the opening time or the closing time:k
tf
=
0.2565
N
·
m.
(3)
The preload torque applied to the valve at the default position is0.3829
N
·
m.Figure 9
Plate angle versus load torque.


### 3.6. Viscous Coefficient Identification

In this experiment, we take off the middle gear; then the motor is separated from the rest of the assembly. When a DC voltage is applied to the motor, the motor can run freely which is given by(16)
E
a
=
R
a
i
a
+
k
b
θ
˙
m
,
J
m
θ
¨
m
=
k
t
i
a
-
k
m
θ
˙
m
.When the motor reaches steady velocity, (16) can be written in form:
(17)
E
a
-
R
a
i
a
k
t
i
a
=
k
b
k
m
.Motor voltage along with the armature current is plotted using Matlab cftool shown in Figure10.Figure 10
Data fit curve.The slope of the line gives the spring constant which readsk
m
=
1.5
×
1
0
-
5
N
·
m
/
rad
/
s. The viscous coefficient of the plate is small and can be ignored in this section. We will identify it in the optimization section.


### 3.7.  Moment of Inertia Identification

Whenθ
˙
>
0 and θ
>
θ
0, the system can be simplified to
(18)
J
θ
¨
+
n
2
k
m
+
k
f
θ
˙
+
k
tf
+
k
sp
θ
-
θ
0
+
k
pre
=
n
k
t
i
a
.Figure11 shows the captured load torque and angular acceleration, calculated by taking derivative of angular velocity, and angular velocity is the derivative of the position. The moment of inertia at the steady time of 13.8 s is 0.0225
kg
·
m
2.Figure 11
Angular acceleration versus load torque.


### 3.8. Optimization

Simulink Design Optimization is a toolbox provided by Matlab, which can estimate the unknown parameters in a model. Due to the approximation of the previous experiments, there are some errors to the model parameters. Therefore, we optimized the previous parameters using parameter estimate tool. In the process of optimization, the initial values are set to the values identified above, and the range for each parameter is defined. The algorithm used on optimization is nonlinear least square method. The optimized parameters are listed in Table1.Table 1
Parameters value.
Parameters
Identified value
Initial
Optimized
J (kg·m2)
0.0225
0.00297
k
b (V·s/rad)
0.0184
0.01681
k
cl
⁡ (Nm/rad)
7.3653
3.6221
k
f (Nm·s/rad)
—
0.00134
k
m (Nm·s/rad)
1.5
e
-
05
1.42130
e
-
05
k
op (Nm/rad)
0.4387
0.21249
k
pre (Nm)
0.3829
0.5306
k
t (Nm·A−1)
0.0184
0.02505
k
tf (Nm)
0.2565
0.37039
R
a (Ω)
2
1.9306
θ
0 (rad)
0.3053
0.302
n
20.68
20.68


### 3.9. Parameters Verification

In this section, we will compare the simulation results with the real measurement. Figure12 shows the response of the real system versus model when subjected to different control inputs. Obviously, the simulations give a good fit to the real outputs.Parameters verification subjected to sinusoidal and step input.
(a)
Subjected to sinusoidal input
(b)
Subjected to a step input


## 3.1.  xPC Target

The real-time simulation of rapid control prototyping (RCP) is that the real plant is controlled by simulating model. For our RCP, we have used xPC Target system (Figure4). Various tests are conducted to estimate the parameters on the platform.Figure 4
xPC Target platform.


## 3.2. Characteristic Curve

The plant has two potentiometers supplied by 5-volt DC. The calibration curves of Throttle Position Sensor (TPS) 1 and TPS2 are shown in Figure5. The angle of the valve is between 10.5° and 101.5° due to the physical limit. When there is no voltage applied to the driven motor, the output voltage of TPS1 is 0.83 V, so the default position θ
0 is 17.5°.Figure 5
Characteristic curve of angle sensors.


## 3.3. Identification ofR
a

There is nearly no rotation when a slowly changing voltage is applied to the motor at the almost closed position [19]. Since the motor current changes slowly, the voltage across the inductance can be neglected. Besides, the inner resistance of the storage can be ignored. Therefore, the armature circuit equations (1) and (2) can be described by
(13)
E
a
=
R
a
i
a
.Figure6 displays the captured curves of the motor voltage against the motor current. Due to the good symmetry of the H-bridge, the direction of the current and voltage is neglected. The inverse of the slope represents the armature resistance which reads R
a
=
2.0
Ω.Figure 6
Motor current versus voltage.


## 3.4. Identification ofk
b

The circuit equations (1) and (2) can be rewritten as
(14)
V
b
=
n
k
b
θ
˙
=
E
a
-
R
a
i
a
,
where n
=
20.68 is the total gear ratio of the two set reduction gear. The back-EMF voltage is generated when the valve shaft is rotating. Figure 7 shows the captured back-EMF voltage and the plate’s angular velocity.Figure 7
Velocity versus back-EMF voltage.We get the back electromotive force constantk
b
=
0.0184 V/(rad/s) at the time when the plate’s angular velocity is steady. Usually, k
b and k
t are almost the same for one motor.


## 3.5.  Static Load Test

Figure8 shows the typical system response when inputting sinusoidal signal. If the frequency is low, the angular velocity and angular acceleration can be ignored. Therefore, the dynamic equation of the system can be described by
(15)
k
sp
θ
-
θ
0
+
k
pre
sgn
⁡
θ
-
θ
0
+
k
tf
sgn
⁡
θ
˙
=
n
k
t
i
a
=
T
.Figure 8
Static loads test.The response plotted in load torque versus angle characteristic is shown in Figure9. This experiment identifies the parameters: k
sp, k
pre, and k
tf.(1)
The inverse slope of the inserted straight lines represented the spring stiffness which readsk
op
=
0.4387
N
·
m
/
rad, k
cl
⁡
=
7.3693
N
·
m
/
rad.
(2)
The friction torque is half of the deviation between two inserted straight lines during the opening time or the closing time:k
tf
=
0.2565
N
·
m.
(3)
The preload torque applied to the valve at the default position is0.3829
N
·
m.Figure 9
Plate angle versus load torque.


## 3.6. Viscous Coefficient Identification

In this experiment, we take off the middle gear; then the motor is separated from the rest of the assembly. When a DC voltage is applied to the motor, the motor can run freely which is given by(16)
E
a
=
R
a
i
a
+
k
b
θ
˙
m
,
J
m
θ
¨
m
=
k
t
i
a
-
k
m
θ
˙
m
.When the motor reaches steady velocity, (16) can be written in form:
(17)
E
a
-
R
a
i
a
k
t
i
a
=
k
b
k
m
.Motor voltage along with the armature current is plotted using Matlab cftool shown in Figure10.Figure 10
Data fit curve.The slope of the line gives the spring constant which readsk
m
=
1.5
×
1
0
-
5
N
·
m
/
rad
/
s. The viscous coefficient of the plate is small and can be ignored in this section. We will identify it in the optimization section.


## 3.7.  Moment of Inertia Identification

Whenθ
˙
>
0 and θ
>
θ
0, the system can be simplified to
(18)
J
θ
¨
+
n
2
k
m
+
k
f
θ
˙
+
k
tf
+
k
sp
θ
-
θ
0
+
k
pre
=
n
k
t
i
a
.Figure11 shows the captured load torque and angular acceleration, calculated by taking derivative of angular velocity, and angular velocity is the derivative of the position. The moment of inertia at the steady time of 13.8 s is 0.0225
kg
·
m
2.Figure 11
Angular acceleration versus load torque.


## 3.8. Optimization

Simulink Design Optimization is a toolbox provided by Matlab, which can estimate the unknown parameters in a model. Due to the approximation of the previous experiments, there are some errors to the model parameters. Therefore, we optimized the previous parameters using parameter estimate tool. In the process of optimization, the initial values are set to the values identified above, and the range for each parameter is defined. The algorithm used on optimization is nonlinear least square method. The optimized parameters are listed in Table1.Table 1
Parameters value.
Parameters
Identified value
Initial
Optimized
J (kg·m2)
0.0225
0.00297
k
b (V·s/rad)
0.0184
0.01681
k
cl
⁡ (Nm/rad)
7.3653
3.6221
k
f (Nm·s/rad)
—
0.00134
k
m (Nm·s/rad)
1.5
e
-
05
1.42130
e
-
05
k
op (Nm/rad)
0.4387
0.21249
k
pre (Nm)
0.3829
0.5306
k
t (Nm·A−1)
0.0184
0.02505
k
tf (Nm)
0.2565
0.37039
R
a (Ω)
2
1.9306
θ
0 (rad)
0.3053
0.302
n
20.68
20.68


## 3.9. Parameters Verification

In this section, we will compare the simulation results with the real measurement. Figure12 shows the response of the real system versus model when subjected to different control inputs. Obviously, the simulations give a good fit to the real outputs.Parameters verification subjected to sinusoidal and step input.
(a)
Subjected to sinusoidal input
(b)
Subjected to a step input


## 4. Controller Design

As mentioned above, the BPV system is featured with strong nonlinearity and two order dynamics. Thus, the simple application of standard PID controller may not work very effectively. ADRC is a new control algorithm based on the active disturbance rejection concept, which can cope with the highly nonlinear dynamics. With ADRC, the unknown dynamics and disturbances can be actively estimated and compensated, and it makes the feedback controller more robust. In fact, system’s equation (7) can be abstracted as the following second-order nonlinear system with external disturbances:
(19)
x
˙
1
=
x
2
,
x
˙
2
=
f
x
1
,
x
2
,
t
+
w
t
+
b
0
u
,
y
=
x
1
,
where x
1, x
2 are system states, f is a smooth function, w
t is the external disturbance, b
0 is the system parameter, and u is the control input signal. The key of ADRC is to design an extended state observer (ESO) to estimate not only states but also the “total disturbance” which contains internal uncertain dynamics and external disturbance.ADRC consists of three parts: transient profile generator, nonlinear state error feedback control law, and extended state observer. Taking second-order controlled object, for example, the structure of active disturbance rejection controller is shown in Figure13 [20].Figure 13
The structure of ADRC, wherev
1 is the desired trajectory of the input and v
2 is its derivative. z
1, z
2 are the estimate of the states x
1 and x
2, and z
3 is the estimate of the system disturbance. u is the static or dynamic feedback control law about e
1 and e
2, where e
1
=
v
1
-
z
1, e
2
=
v
2
-
z
2. Since the BPV is a second-order system, therefore, the controller is designed based on the two order control algorithms.


### 4.1. Transient Profile Generator

Reasonable arrangement of the transition process makes the set-value input to the control algorithm not mutate but changes according to a certain expectation. For example, in the classical control theory, predicting filter is often adopted to process input signal. Transient profile generator is a second-order system which can arrange smooth transition process according to the object ability and control necessary. Meanwhile, it collects the differential signals of each order. It is an effective method to solve the contradiction between quickness and overshoot. A discrete-time solution for a discrete double integral plant is described by(20)
v
1
k
+
1
=
v
1
k
+
h
*
v
2
k
,
v
2
k
+
1
=
v
2
k
+
h
*
fhan
e
1
,
v
2
,
r
0
,
h
0
,
T
0
=
2
A
r
0
,
where A is the step magnitude and T
0 is the rising time of the input signal. h is the sampling period and r
0 and h
0 (h
0
=
k
h, k
∈
Z
+ or h
0
=
h) are controllers parameters, which are adjusted accordingly as filter coefficients. See the appendix for undefined function fhan
(
e
1
,
v
2
,
r
0
,
h
0
) (similarly hereinafter).


### 4.2.  Extended State Observer

Systems are working under different kinds of disturbances, among which the ones that have effects on the output signal are the most important. So we can separate them from the output signal by defining a new state. This is done by ESO. The ESO provides the estimate of the unmeasured system’s state and the real-time action of the unknown disturbances and then compensates them. The ESO can improve the performance of the system, which is in the form of(21)
e
1
=
v
1
-
z
1
,
f
1
=
fal
e
1
,
0.5
,
δ
,
f
2
=
fal
e
1
,
0.25
,
δ
,
z
1
k
+
1
=
z
1
k
+
h
z
2
k
-
β
01
e
1
,
z
2
k
+
1
=
z
2
k
+
h
z
3
k
-
β
02
f
1
+
b
^
0
u
,
z
3
k
+
1
=
z
3
k
+
h
-
β
03
f
2
,
where z
1, z
2 are the estimate of the states x
1 and x
2 and z
3 is the estimate of the system disturbances. β
01, β
02, and β
03 are observer gains and b
^
0 represents the nominal value of the system parameter b
0. A great number of simulations show that [20] if δ
=
h, β
01
=
1
/
h then β
02 and β
03 are functions of h.


### 4.3.  Feedback Control

In the feedback control law,u
0 is a function of static or dynamic state errors e
1 and e
2:
(22)
e
1
=
v
1
-
z
1
,
e
2
=
v
2
-
z
2
,
u
0
=
k
e
1
,
e
2
.Disturbance compensation is the following:(23)
u
=
u
0
-
z
3
b
0
.


### 4.4. Parameters Tuning

According to the control targets made by BOSCH shown in Table2 [21], t
r
<
0.1
s, A
=
(
θ
op
-
θ
cl
⁡
)
-
θ
0
*
92
%, then 2
A
/
r
0
<
0.1, so r
0 can be a value of 800. In addition, h
=
t
s
/
100
=
0.001
(
s
), h
0
=
h.Table 2
The performance of step response.
Performance indicators
Position angle
4–96%
96–4%
Rise time (ms)
<100
<100
Regulation time (ms)
<300
<300
Overshoot permitted (%)
<0.6
<3.5Curves in [20] obtained from simulation experiments show that if δ
=
h, β
01
=
1
/
h, then β
02 and β
03 are functions of h. β
01
≈
1
/
1.6
h
1.5, β
02
≈
1
/
8.6
h
2.2. So we can set β
02
=
16000, β
03
=
500000. Besides, the nominal value of the system b
^
0
=
b
0
=
n
k
t
/
J
R
a
=
90.The PID control law is selected as the state error feedback control law for our BPV system, which is described by(24)
u
0
=
K
p
e
1
+
K
i
∫
0
t
e
1
d
t
+
K
d
e
2
.Based on the PID tuning method, the control parameters arek
p
=
20, k
i
=
15, and k
d
=
10.


### 4.5. Simulations and Analyses

The controller described in previous sections has been tested by simulation model and the real BPV system with the identified parameters estimated in Section3. There are two physical limitations in our system, which have been discussed in the previous sections. Then, the maximum opening is 101.5° and the limp home angle is 17.5°. As the performance indicators shown in Table 2, the input signals are designed as firstly, θ
(
t
) steps down from 98° (96%) to 20° (4%), staying unchanged for 1 second, and then steps to 98° (96%). We also design a tracking control input θ
(
t
)
=
59
+
39
sin
π
t.Figures14 and 15 show the simulation results implemented in Simulink. The step responses from 4% to 96% and from 96% to 4% have nearly the same results. Both the rise time and the regulation time of the responses are less than 100 ms, and with no overshoot, which are fit for the technical targets made by BOSCH. From the simulation result of the sine test, one can see that the tracking error is less than 0.2° via ADRC. Besides, the estimate disturbance is in good agreement with the real disturbance for both inputs.Subjected to a step input (96%–4%–96%).
(a)
Step response
(b)
Disturbance estimateSubjected to a sinusoidal input.
(a)
Dynamic errors
(b)
Disturbance estimateThe real-time ADRC is implemented with xPC Target shown in Figure16, where k
p
=
20, k
i
=
5, and k
d
=
0.35. The results show that the controller can satisfy the preconcert control request.xPC Target test.
(a)
Subjected to a step input (4%–96%–4%)
(b)
Subjected to a sinusoidal inputThe test results of the performance are shown in Table3.Table 3
The performance of step response.
Test type
Indicator
Position angle
Rise time (ms)
Regulation time (ms)
Overshoot (%)
Model simulation
4–96%
<70
<90
<0.2
96–4%
<80
<100
<0.5
Real-time test
4–96%
<80
<130
<0.1
96–4%
<80
<150
<0.25


## 4.1. Transient Profile Generator

Reasonable arrangement of the transition process makes the set-value input to the control algorithm not mutate but changes according to a certain expectation. For example, in the classical control theory, predicting filter is often adopted to process input signal. Transient profile generator is a second-order system which can arrange smooth transition process according to the object ability and control necessary. Meanwhile, it collects the differential signals of each order. It is an effective method to solve the contradiction between quickness and overshoot. A discrete-time solution for a discrete double integral plant is described by(20)
v
1
k
+
1
=
v
1
k
+
h
*
v
2
k
,
v
2
k
+
1
=
v
2
k
+
h
*
fhan
e
1
,
v
2
,
r
0
,
h
0
,
T
0
=
2
A
r
0
,
where A is the step magnitude and T
0 is the rising time of the input signal. h is the sampling period and r
0 and h
0 (h
0
=
k
h, k
∈
Z
+ or h
0
=
h) are controllers parameters, which are adjusted accordingly as filter coefficients. See the appendix for undefined function fhan
(
e
1
,
v
2
,
r
0
,
h
0
) (similarly hereinafter).


## 4.2.  Extended State Observer

Systems are working under different kinds of disturbances, among which the ones that have effects on the output signal are the most important. So we can separate them from the output signal by defining a new state. This is done by ESO. The ESO provides the estimate of the unmeasured system’s state and the real-time action of the unknown disturbances and then compensates them. The ESO can improve the performance of the system, which is in the form of(21)
e
1
=
v
1
-
z
1
,
f
1
=
fal
e
1
,
0.5
,
δ
,
f
2
=
fal
e
1
,
0.25
,
δ
,
z
1
k
+
1
=
z
1
k
+
h
z
2
k
-
β
01
e
1
,
z
2
k
+
1
=
z
2
k
+
h
z
3
k
-
β
02
f
1
+
b
^
0
u
,
z
3
k
+
1
=
z
3
k
+
h
-
β
03
f
2
,
where z
1, z
2 are the estimate of the states x
1 and x
2 and z
3 is the estimate of the system disturbances. β
01, β
02, and β
03 are observer gains and b
^
0 represents the nominal value of the system parameter b
0. A great number of simulations show that [20] if δ
=
h, β
01
=
1
/
h then β
02 and β
03 are functions of h.


## 4.3.  Feedback Control

In the feedback control law,u
0 is a function of static or dynamic state errors e
1 and e
2:
(22)
e
1
=
v
1
-
z
1
,
e
2
=
v
2
-
z
2
,
u
0
=
k
e
1
,
e
2
.Disturbance compensation is the following:(23)
u
=
u
0
-
z
3
b
0
.


## 4.4. Parameters Tuning

According to the control targets made by BOSCH shown in Table2 [21], t
r
<
0.1
s, A
=
(
θ
op
-
θ
cl
⁡
)
-
θ
0
*
92
%, then 2
A
/
r
0
<
0.1, so r
0 can be a value of 800. In addition, h
=
t
s
/
100
=
0.001
(
s
), h
0
=
h.Table 2
The performance of step response.
Performance indicators
Position angle
4–96%
96–4%
Rise time (ms)
<100
<100
Regulation time (ms)
<300
<300
Overshoot permitted (%)
<0.6
<3.5Curves in [20] obtained from simulation experiments show that if δ
=
h, β
01
=
1
/
h, then β
02 and β
03 are functions of h. β
01
≈
1
/
1.6
h
1.5, β
02
≈
1
/
8.6
h
2.2. So we can set β
02
=
16000, β
03
=
500000. Besides, the nominal value of the system b
^
0
=
b
0
=
n
k
t
/
J
R
a
=
90.The PID control law is selected as the state error feedback control law for our BPV system, which is described by(24)
u
0
=
K
p
e
1
+
K
i
∫
0
t
e
1
d
t
+
K
d
e
2
.Based on the PID tuning method, the control parameters arek
p
=
20, k
i
=
15, and k
d
=
10.


## 4.5. Simulations and Analyses

The controller described in previous sections has been tested by simulation model and the real BPV system with the identified parameters estimated in Section3. There are two physical limitations in our system, which have been discussed in the previous sections. Then, the maximum opening is 101.5° and the limp home angle is 17.5°. As the performance indicators shown in Table 2, the input signals are designed as firstly, θ
(
t
) steps down from 98° (96%) to 20° (4%), staying unchanged for 1 second, and then steps to 98° (96%). We also design a tracking control input θ
(
t
)
=
59
+
39
sin
π
t.Figures14 and 15 show the simulation results implemented in Simulink. The step responses from 4% to 96% and from 96% to 4% have nearly the same results. Both the rise time and the regulation time of the responses are less than 100 ms, and with no overshoot, which are fit for the technical targets made by BOSCH. From the simulation result of the sine test, one can see that the tracking error is less than 0.2° via ADRC. Besides, the estimate disturbance is in good agreement with the real disturbance for both inputs.Subjected to a step input (96%–4%–96%).
(a)
Step response
(b)
Disturbance estimateSubjected to a sinusoidal input.
(a)
Dynamic errors
(b)
Disturbance estimateThe real-time ADRC is implemented with xPC Target shown in Figure16, where k
p
=
20, k
i
=
5, and k
d
=
0.35. The results show that the controller can satisfy the preconcert control request.xPC Target test.
(a)
Subjected to a step input (4%–96%–4%)
(b)
Subjected to a sinusoidal inputThe test results of the performance are shown in Table3.Table 3
The performance of step response.
Test type
Indicator
Position angle
Rise time (ms)
Regulation time (ms)
Overshoot (%)
Model simulation
4–96%
<70
<90
<0.2
96–4%
<80
<100
<0.5
Real-time test
4–96%
<80
<130
<0.1
96–4%
<80
<150
<0.25


## 5. Conclusion

In this paper, the mathematical modeling of the BPV has been further explored and identification algorithms of proposed model parameters have been designed based on series of measurements of the real system responses using xPC Target platform. Due to some unique features in back-pressure valve model structures, parameters identification encounters some challenging issues. Different identification algorithms are introduced to resolve these issues. After that, the ADRC algorithm is introduced. The key component is reference profile generation, which can avoid set point jump, so that the output of the plant can follow reasonably. The second important part of the controller is the extended state observe. Both the model error and external disturbance are regarded as additional external disturbance and they are estimated and compensated by the extended state observers in real time.Finally, the controller has been implemented in the simulation model and the actual system. The simulations show that the proposed controller can not only achieve a good dynamic and steady performance but also have strong robustness against parameters’ disturbance and external disturbance. Extensive researches will be done to combine the BPV model with the PEMFC system model to obtain optimized opening at a given operation.



---

*Source: 246015-2014-11-16.xml*