Implement generic hydrogen fuel cell stack model

Simscape / Electrical / Specialized Power Systems / Electric Drives / Extra Sources

The Fuel Cell Stack block implements a generic model parameterized to represent most popular types of fuel cell stacks fed with hydrogen and air.

The block represents two versions of the stack model: a simplified model and a detailed
model. You can switch between the two models by selecting the level in the mask under
**Model detail level** in the block dialog box.

This model is based on the equivalent circuit of a fuel cell stack shown below:

The simplified model represents a particular fuel cell stack operating at nominal conditions of temperature and pressure. The parameters of the equivalent circuit can be modified based on the polarization curve obtained from the manufacturer datasheet. You just have to input in the mask the value of the voltage at 0 and 1 A, the nominal and the maximum operating points, for the parameters to be calculated. A diode is used to prevent the flow of negative current into the stack. A typical polarization curve consists of three regions:

The first region represents the activation voltage drop due to the slowness of the chemical reactions taking place at electrode surfaces. Depending on the temperature and operating pressure, type of electrode, and catalyst used, this region is more or less wide. The second region represents the resistive losses due to the internal resistance of the fuel cell stack. Finally, the third region represents the mass transport losses resulting from the change in concentration of reactants as the fuel is used.

The detailed model represents a particular fuel cell stack when the parameters such as
pressures, temperature, compositions and flow rates of fuel and air vary. You can select which
parameters to vary on the **Signal variation** pane on the block dialog box.
These variations affect the open circuit voltage
(*E*_{oc}), the exchange current
(*i*_{0}), and the Tafel slope
(*A*). *E*_{oc},
*i*_{0} and *A*
are modified as follows:

$$\begin{array}{c}{E}_{oc}={K}_{c}{E}_{n}\\ {i}_{0}=\frac{zFk\left({P}_{{H}_{2}}+{P}_{{O}_{2}}\right)\Delta v}{Rh}{e}^{\frac{-\Delta G}{RT}}\\ A=\frac{RT}{z\alpha F},\end{array}$$

where:

*R* = 8.3145 J/(mol K)

*F* = 96485 A s/mol

*z* = Number of moving electrons

*E*_{n} = Nernst voltage, which
is the thermodynamics voltage of the cells and depends on the temperatures and partial pressures
of reactants and products inside the stack (V)

α = Charge transfer coefficient, which depends on the type of electrodes and catalysts used

*P*_{H2} = Partial pressure of
hydrogen inside the stack (Pa)

*P*_{O2} = Partial pressure of
oxygen inside the stack (Pa)

*k* = Boltzmann's constant = 1.38 ×
10^{–23} J/K

*h* = Planck's constant = 6.626 × 10^{–34}
J s

*Δv* = Activation barrier volume factor (m^{3}).
The size of activation barrier (*ΔG*) is computed assuming
*Δv* = 1 m^{3}.

*ΔG* = Size of the activation barrier which depends on the type of
electrode and catalyst used (J/mol)

*T* = Temperature of operation (K)

*K*_{c} = Voltage constant at
nominal condition of operation

The equivalent circuit is the same as for the simplified model, except that the parameters
*E*_{oc},
*i*_{0} and *Α*
have to be updated on-line as shown below:

The rates of conversion (utilizations) of hydrogen
(*U*_{fH2}) and oxygen
(*U*_{fO2}) are determined in Block
A as follows:

$$\{\begin{array}{c}{U}_{f{H}_{2}}=\frac{{n}_{{H}_{2}}^{r}}{{n}_{{H}_{2}}^{\text{in}}}=\frac{60000RTN{i}_{fc}}{zF{P}_{\text{fuel}}{V}_{lpm(\text{fuel})}x\%}\\ {U}_{f{O}_{2}}=\frac{{n}_{{O}_{2}}^{r}}{{n}_{{O}_{2}}^{\text{in}}}=\frac{60000RTN{i}_{fc}}{2zF{P}_{\text{air}}{V}_{lpm(\text{air})}y\%}\end{array}$$

where

*P*_{fuel} = Absolute supply
pressure of fuel (atm)

*P*_{air} = Absolute supply
pressure of air (atm)

*V*_{lpm(fuel)} = Fuel flow
rate (l/min)

*V*_{lpm(air)} = Air flow rate
(l/min)

*x* = Percentage of hydrogen in the fuel (%)

*y* = Percentage of oxygen in the oxidant (%)

*N* = Number of cells

The 60000 constant comes from the conversion from the liter/min flow rate used in the model to m3/s (1 liter/min = 1/60000 m3/s).

The partial pressures and the Nernst voltage are determined in Block B as follows:

$$\{\begin{array}{c}{P}_{{H}_{2}}=\left(1-{U}_{{f}_{{H}_{2}}}\right)x\%{P}_{\text{fuel}}\\ {P}_{{H}_{2}O}=\left(w+2y\%{U}_{{f}_{{O}_{2}}}\right){P}_{\text{air}}\\ {P}_{{O}_{2}}=\left(1-{U}_{{f}_{{O}_{2}}}\right)y\%{P}_{\text{air}}\end{array}$$

and

$${E}_{n}=\{\begin{array}{cc}1.229+(T-298)\frac{-44.43}{zF}+\frac{RT}{zF}\mathrm{ln}\left({P}_{{H}_{2}}{P}_{{O}_{2}}^{1/2}\right)& \text{when}T\le {100}^{\circ}C\\ 1.229+(T-298)\frac{-44.43}{zF}+\frac{RT}{zF}\mathrm{ln}\left(\frac{{P}_{{H}_{2}}{P}_{{O}_{2}}^{1/2}}{{P}_{{H}_{2}O}}\right)& \text{when}T{100}^{\circ}C\end{array}$$

where

*P*_{H2O} = Partial pressure of
water vapor inside the stack (atm)

*w* = Percentage of water vapor in the oxidant (%)

From the partial pressures of gases and the Nernst voltage, the new values of the open
circuit voltage (*E*_{oc}) and the
exchange current (*i*_{0}) can be
calculated.

Block C calculates the new value of the Tafel slope (A).

The parameters α, *ΔG*, and
*K*_{c} are calculated based on
the polarization curve at nominal conditions of operation along with some additional parameters,
such as the low heating value (LHV) efficiency of the stack, composition of fuel and air, supply
pressures and temperatures. They can be easily obtained from the manufacturer datasheet.

The nominal rates of conversion of gases are calculated as follows:

$$\begin{array}{c}{U}_{{f}_{{H}_{2}}}=\frac{{\eta}_{\text{nom}}\Delta {h}^{0}\left({H}_{2}O(\text{gas})\right)N}{zF{V}_{\text{nom}}}\\ {U}_{{f}_{{O}_{2}}}=\frac{60000R{T}_{\text{nom}}N{I}_{\text{nom}}}{2zF{P}_{{\text{air}}_{\text{nom}}}{V}_{lpm{(\text{air})}_{\text{nom}}}\cdot 0.21}\end{array}$$

where:

*η*_{nom }= Nominal LHV
efficiency of the stack (%).

*Δh*^{0}*(H*_{2}*O(gas))*
= 241.83 × 10^{3} J/mol.

*V*_{nom} = Nominal voltage
(V).

*I*_{nom} = Nominal current
(A).

*V*_{lpm(air)nom} = Nominal air
flow rate (l/min).

*P*_{airnom} = Nominal absolute
air supply pressure (Pa).

*T*_{nom} = Nominal operating
temperature (K).

From these rates of conversion, the nominal partial pressures of gases and the Nernst
voltage can be derived. With
*E*_{oc},
*i*_{0} and *Α*
known and assuming that the stack operates at constant rates of conversion or utilizations at
nominal condition, α, *ΔG*, and
*K*_{c} can be determined.

If there is no fuel or air at the stack input, it is assumed that the stack is operating at a fixed rate of conversion of gases (nominal rate of conversion), that is, the supply of gases is adjusted according to the current so that they are always supplied with just a bit more than needed by the stack at any load.

The maximum current the stack can deliver is limited by the maximum flow rates of fuel and air that can be reached. Beyond that maximum current, the voltage output by the stack decreases abruptly as more current is drawn.

The dynamics of the fuel cell are represented if you specify the response time and the
parameters for flow dynamics (peak utilization and corresponding voltage undershoot) on the
**Fuel Cell Dynamics** pane on the dialog box.

The response time (*T*_{d}) @
95% is used to model the "charge double layer" phenomenon due to the build-up of charges at
electrode/electrolyte interface. This affects only the activation voltage
(*NA*ln*(i*_{fc}*/i*_{0}*)*)
as shown on the equivalent circuits.

The peak utilization
*(Uf*_{O2(peak)}) and the
corresponding voltage undershoot
(*V*_{u}) are used to model the
effect of oxygen depletion (due to the air compressor delay) on the cell output voltage. The
Nernst voltage is modified due to this effect as follows:

$${E}_{n}=\{\begin{array}{ll}{E}_{n}-K\left({U}_{{f}_{{O}_{2}}}-{U}_{{f}_{{O}_{2}}(\text{nom})}\right)\hfill & {U}_{{f}_{{O}_{2}}}>{U}_{{f}_{{O}_{2}}(\text{nom})}\hfill \\ {E}_{n}\hfill & {U}_{{f}_{{O}_{2}}}\le {U}_{{f}_{{O}_{2}}(\text{nom})}\hfill \end{array}$$

where

*K =* voltage undershoot constant

*U*_{fO2(nom)} = nominal oxygen
utilization

*K* is determined as follows:

$$K=\frac{{V}_{u}}{{K}_{c}\left({U}_{{f}_{{O}_{2}}(\text{peak})}-{U}_{{f}_{{O}_{2}}(\text{nom})}\right)}.$$

Current step and interrupt tests must be made on a real stack to represent with accuracy
its dynamics. The figure below shows the stack response from these tests and the required
parameters (*T*_{d},
*Uf*_{O2(peak)} and
*V*_{u}).

The response time (*T*_{d})
depends on the fuel cell stack itself and is usually given on the datasheet. The parameters for
flow dynamics (*Uf*_{O2(peak)} and
*V*_{u}) depend on the dynamics of
external equipment (compressor, regulator, and loads) and they are not provided by manufacturers
as their values vary with the user application. For simulation, assume values of
*Uf*_{O2(peak)} between 60% to 70%
and *V*_{u} between 2-5% of the
stack nominal voltage.

**Preset model**Provides a set of predetermined polarization curves and parameters for particular fuel cell stacks found on the market:

`No (User-Defined)`

(default)`PEMFC - 1.26 kW - 24 Vdc`

`PEMFC - 6 kW - 45 Vdc`

`PEMFC - 50 kW - 625 Vdc`

`AFC - 2.4 kW - 48 Vdc`

`SOFC - 3 kW - 100 Vdc`

`SOFC - 25 kW - 630 Vdc`

Select one of these preset models to load the corresponding parameters in the entries of the dialog box. Select

`No (User-Defined)`

if you do not want to use a preset model.**Model detail level**Provide access to the two versions of the model:

`Simplified`

`Detailed`

(default)

When a simplified model is used, there is no variable under the

**signal variation**tab.**Voltage at 0 A and 1 A**The voltage at 0 A and 1 A of the stack (Volts). Assuming nominal and constant gases utilizations. Default is

`[65 63]`

.**Nominal operating point**The rated current (Ampere) and rated voltage (Volts) of the stack. Assuming nominal and constant gases utilizations. Default is

`[133.3 45]`

.**Maximum operating point**The current (Ampere) and voltage (Volts) of the stack at maximum power. Assuming nominal and constant gases utilizations. Default is

`[225 37]`

.**Number of cells**The number of cells in series in the stack. This parameter is available only for a detailed model. Default is

`65`

.**Nominal stack efficiency**The rated efficiency of the stack relative to the low heating value (LHV) of water. This parameter is available only for a detailed model. Default is

`55`

.**Operating temperature**The nominal temperature of operation in degrees Celsius. This parameter is available only for a detailed model. Default is

`65`

.**Nominal air flow rate**The rated air flow rate (l/min). This parameter is available only for a detailed model. Default is

`300`

.**Nominal supply pressure**Rated supply pressure (absolute) of fuel and air in bars. This parameter is available only for a detailed model. Default is

`[1.5 1]`

.**Nominal composition**The rated percentage of hydrogen (x) in the fuel, oxygen (y) and water (w) in the oxidant. This parameter is available only for a detailed model. Default is

`[99.95 21 1]`

.**Plot V-I characteristic**Plots a figure containing two graphs. The first graph represents the stack voltage (Volts) vs current (A) and the second graph represents the stack power (kW) vs current (A). This button is available only for a detailed model.

**View cell parameters**Presents the overall parameters of the stack. This button is available only for a detailed model. The dialog box is shown below.

The **Signal Variation** tab is available only if **Model detail
level** is set to `Detailed`

. It provides a list of
parameters that can be varied. Select a check box for a variable to input a corresponding
signal to the block. The following signals can be input to the block:

**Fuel composition**Percentage of hydrogen in the fuel. Default is cleared.

**Oxidant composition**Percentage of oxygen in the oxidant. Default is cleared.

**Fuel flow rate**Fuel flow rate in liter per minutes. Default is cleared.

**Air flow rate**Air flow rate in liter per minutes. Default is cleared.

**System Temperature**Temperature of operation in Kelvin. Default is cleared.

**Fuel supply pressure**Fuel supply pressure in bars. Default is cleared.

**Air supply pressure**Air supply pressure in bars. Default is cleared.

**Specify Fuel Cell Dynamics?**Asks whether you want to specify the fuel cell dynamics. Select the check box to enter the fuel cell response time in seconds. Default is cleared.

**Fuel Cell response time (sec)**Enter the response time of the cell (at 95% of the final value). Default is

`1`

. This parameter becomes available only when the**Specify Fuel Cell Dynamics?**check box is selected.**Peak O2 utilization (%)**Enter the peak oxygen utilization at nominal condition of operation. Default is

`80`

. This parameter is available only when the**Specify Fuel Cell Dynamics?**check box is selected and the**Air flow rate**check box is selected in the**Signal Variation**tab.**Voltage undershoot (V) @ peak O2 utilization**Enter the voltage undershoot (Volts) at peak oxygen utilization at nominal condition of operation. Default is

`10`

. This parameter is available only when the**Specify Fuel Cell Dynamics?**check box is selected and the**Air flow rate**check box is selected in the**Signal Variation**tab.

Here is the procedure to extract parameters from fuel cell stack manufacturer's data sheet. For this example, the NetStack PS6 data sheet from NetStack is used:

The rated power of the stack is 6 kW and the nominal voltage is 45 V. The following detailed parameters are deduced from the datasheet.

Voltage at 0 A and 1 A [E

_{oc},V_{1}] = [65, 63]Nominal operating point [I

_{nom}, V_{nom}]=[133.3, 45]Maximum operating point [I

_{end},V_{end}]=[225, 37]Nominal stack efficiency (

*η*_{nom})= 55%Operating temperature = 65

^{0}CNominal supply pressure [H2, Air]=[1.5 1]

If the pressure given is relative to the atmospheric pressure, add 1 bar to get the absolute pressure.

Nominal composition (%)[H2, O2, H2O(Air)]=[99.999, 21, 1]

If air is used as oxidant, assume 21% of O

_{2}and 1% of H_{2}O in case their percentages are not specified.Number of cells

If not specified, estimate it from the formulae below:

$$N=\frac{2\cdot 96485\cdot {V}_{\text{nom}}}{241.83\cdot {10}^{3}\cdot {\eta}_{\text{nom}}}.$$

In this case,

$$N=\frac{2\cdot 96485\cdot 45}{241.83\cdot {10}^{3}\cdot 0.55}=65.28=65\text{cells}\text{.}$$

Nominal air flow rate

If the maximum air flow rate is given, the nominal flow rate can be calculated assuming a constant oxygen utilization at any load. The current drawn by the cell is linearly dependent on air flow rate and the nominal flow rate is given by:

$${V}_{lpm{(\text{air})}_{\text{nom}}}=\frac{{I}_{\text{nom}}\cdot {V}_{lpm{(\text{air})}_{\mathrm{max}}}}{{I}_{\text{end}}}.$$

In this case,

$${V}_{lpm{(\text{air})}_{\text{nom}}}=\frac{133.3\cdot 500}{225}=297\text{liters/min}.$$

In case no information is given, assume the rate of conversion of oxygen to be 50% (as it is usually the case for most fuel cell stacks) and use the formulae below to determine the nominal air flow rate.

$${V}_{lpm{(\text{air})}_{\text{nom}}}=\frac{60000R{T}_{\text{nom}}N{I}_{\text{nom}}}{2zF{P}_{{\text{air}}_{\text{nom}}}\cdot 0.5\cdot 0.21}.$$

Fuel cell response time = 10 s

The parameters [E_{oc}, V_{1}],
[I_{nom}, V_{nom}], and
[I_{end},V_{end}] are approximate and depend on the
precision of the points obtained from the polarization curve. The higher the accuracy of these
parameters, the more closed the simulated stack voltage to the data sheet curve is. A tool,
called ScanIt (from amsterchem) can be used to extract precise values from data sheet
curves.

With the above parameters, the polarization curve of the stack operating at fixed nominal rate of conversion of gases is closed to the datasheet curves as shown below: The blue dotted line shows that the simulated stack voltage and green dotted line shows the simulated stack power.

Above the maximum current, the flow rate of gases entering the stack is maximum and the stack voltage decreases abruptly as more current is drawn.

**m**The Simulink

^{®}output of the block is a vector containing 11 signals. You can demultiplex these signals by using the Bus Selector block provided in the Simulink library.Model detail level

(Signal availability)Signal

Definition

Units

Symbol

Detailed

Simplified

1

Voltage

V

V

_{fc}Yes

Yes

2

Current

I

I

_{fc}Yes

Yes

3

Stack Efficiency

%

η

Yes

No(Set to 0)

4

Stack consumption [Air, Fuel]

slpm

V

_{slpm}Yes

No(Set to 0)

5

Flow Rate [Air, Fuel]

lpm

Fr

_{lpm}Yes

No(Set to 0)

6

Stack consumption [Air, Fuel]

lpm

V

_{lpm}Yes

No(Set to 0)

7

Utilization [Oxygen, Hydrogen]

%

U

_{f}Yes

No(Set to 0)

8

Slope of Tafel curve

A

Yes

No(Set to 0)

9

Exchange current

A

i

_{0}Yes

No(Set to 0)

10

Nernst voltage

V

E

_{n}Yes

No(Set to 0)

11

Open circuit voltage

V

E

_{oc}Yes

No(Set to 0)

The gases are ideal

The stack is fed with hydrogen and air

The stack is equipped with a cooling system which maintains the temperature at the cathode and anode exits stable and equal to the stack temperature

The stack is equipped with a water management system to maintain the humidity inside the cell at appropriate level at any load

The cell voltage drops are due to reaction kinetics and charge transport as most fuel cells do not operate in the mass transport region

Pressure drops across flow channels are negligible

The cell resistance is constant at any condition of operation

Chemical reaction dynamics caused by partial pressure changes of chemical species inside the cell are not considered

The stack output power is limited by the fuel and air flow rates supplied

The effect of temperature and humidity of the membrane on the internal resistance is not considered

The flow of gases or water through the membrane is not considered

The `power_fuel_cell`

example illustrates a 6 kW, 45 volts Proton Exchange Membrane (PEM) Fuel Cell Stack model
feeding a 100Vdc DC/DC converter.

[1] Njoya, S. M., O. Tremblay, and L. -A. Dessaint. *A generic fuel cell
model for the simulation of fuel cell vehicles*. Vehicle Power and Propulsion
Conference, 2009, VPPC ’09, *
IEEE ^{®}
*. Sept. 7–10, 2009, pp. 1722–29.

[2] Motapon, S.N., O. Tremblay, and L. -A. Dessaint.
“Development of a generic fuel cell model: application to a fuel cell vehicle
simulation.” *Int. J. of Power Electronics*. Vol. 4, No. 6, 2012,
pp. 505–22.