?
?
MODELING AND ANALYSIS OF A PEM FUEL CELL SYSTEM  
FOR A QUADRUPED ROBOT 
Except where reference is made to the work of others, the work described in this thesis is 
my own or was done in collaboration with my advisory committee. This thesis does not 
include proprietary or classified information. 
 
_____________________________________ 
Heon Joong Lee 
 
Certificate of Approval:                                                          
 
 
      
Jay Khodadadi       Song-Yul Choe, Chair 
Professor        Associate Professor 
Mechanical Engineering        Mechanical Engineering 
 
      
Sushil Bhavnani      George Flowers 
Professor       Dean of the Graduate School 
Mechanical Engineering     Graduate School 
 
  
?
?
MODELING AND ANALYSIS OF A PEM FUEL CELL SYSTEM  
FOR A QUADRUPED ROBOT 
 
 
Heon Joong Lee 
 
 
A Thesis 
Submitted to 
the Graduate Faculty of 
Auburn University 
in Partial Fulfillment of the 
Requirements for the 
Degree of 
Master of Science 
 
 
Auburn, Alabama 
August 10, 2009 
iii?
?
MODELING AND ANALYSIS OF A PEM FUEL CELL SYSTEM 
FOR A QUADRUPED ROBOT 
 
 
 
Heon Joong Lee 
 
 
 
Permission is granted to Auburn University to make copies of this thesis at its discretion, 
upon the request of individuals or institutions and at their expense. 
The author reserves all publication rights. 
 
 
 
             ________________________ 
         Signature of Author 
 
  
       ________________________ 
             Date of Graduation 
iv?
?
VITA 
Heon Joong Lee, son of Dae Kyu Lee and Kyung Ae You, was born on Jan 10
th
, 
1981 in Seoul, Korea. He earned a Bachelor of Science degree in the Department of 
Mechanical Engineering at Sogang University in Feb, 2007. He has studied at Auburn 
University pursuing a Master of Science degree in the Department of Mechanical 
Engineering since May 2007. 
??
v?
?
THESIS ABSTRACT 
MODELING AND ANALYSIS OF A PEM FUEL CELL SYSTEM  
 FOR A QUADRUPED ROBOT 
 
Heon Joong Lee 
Master of Science, August 10, 2009 
(B.S., Sogang University, 2007) 
135 Typed Pages 
Directed by Song-Yul Choe 
 
The polymer electrolyte membrane fuel cell is a prospective alternative energy 
conversion device that has environment-friendly characteristics. It has not only 
electrochemical reaction but also coupled thermal fluid characteristics. Chemical reaction 
of the reactants, mass transport of products, heat generation and heat transfer from 
catalyst layer to the end plate occur during operation of the fuel cell. In particular, water 
management and heat management play an important role in safe and successful 
operation of the Polymer Electrolyte Membrane (PEM) fuel cell system. Therefore, the 
understanding of the dynamic characteristics of the PEM fuel cell stack is crucial to 
ensure accurate analysis and effective design of the PEM fuel cell system.  
Generally, computational methodology is an essential tool for investigating 
interrelated physical parameters and dynamic behavior. The dynamic model for the PEM
vi?
?
 fuel cell system was developed considering temperature and two-phase effects in the 
Advanced Propulsion Research Laboratory at Auburn University. The model developed 
and compared with experimental results. Moreover, the two cell stack model is used to 
analyze static and dynamic behaviors under different operating conditions that include 
operating temperature, relative humidity of air, stoichiometric ratio, partial pressure of 
reactants in addition to various geometric parameters that include thickness and porosity 
of the gas diffusion layer (GDL) and the thickness of the membrane under the ramp and 
step current loads.  
The analysis of the 20-cell stack is also performed to investigate the temperature 
effect between cells within the stack before a multi cell stack is constructed. The analysis 
consists of static and dynamic behavior for understanding temperature distribution and 
voltage difference due to the temperature effect. Finally, the 100-cell stack is constructed 
and its characteristics are analyzed. However, stack model alone has a limitation for 
estimating real net power from the fuel cell system without consideration of auxiliary 
components such as an air processing system and a thermal circuit. Therefore, each 
component is modeled for the 100-cell PEM fuel cell system. 
 Various comparative studies showed that proper operating conditions are critical 
for improved performance and two phase effects which both affect concentration of 
reactants, thus should be considered for the design and optimization of the stack. 
Moreover, appropriate design and sizing of the air supply system and thermal circuit 
components could reduce parasitic power and increased efficiency of the fuel cell system. 
In the multi cell stack model, the temperature distribution affected not only concentration 
vii?
?
of the reactant but also proton conductivity of the membrane. Consequently, it caused the 
change in concentration overpotential, ohmic overpotential and cell voltage. 
viii?
?
ACKNOWLEDGEMENTS 
I would like to express special gratitude to the advisor and chair of the committee, 
Professor Song-Yul Choe, for his guidance during my time here at Auburn University. 
I would also like to thank the other committee members, Professor Jay Khodadadi 
and Professor Sushil Bhavnani for their valuable advice. Special thanks are extended to 
all of my colleagues, Jinglin He, Kevin Siniard, Alan Dunlavy, Xiao Meng, Ruijian Fu 
and Jong-Woo Ahn on the Advanced Propulsion Laboratorty of Auburn University for 
constructive discussions during research work. 
Lastly, I would like to thank my parents and wife, Dae-Kyu Lee, Kyung-Ae Ryu 
and Sunja Kim. Without their strong motivation and support, I would not be where I am 
today. 
 
 
  
ix?
?
Computer software used: Microsoft word 2007 
                    Matlab/Simulink 7.0.4 
                    Origin 7.5
x?
?
TABLE OF CONTENTS 
LIST OF FIGURES ......................................................................................................... xiii?
LIST OF TABLES .......................................................................................................... xvii?
NOMENCLATURE ...................................................................................................... xviii?
1? INTRODUCTION .................................................................................................. 1?
1.1? THE FUEL CELL ................................................................................................... 1?
1.2? PRINCIPLES OF A FUEL CELL OPERATION .................................................. 3?
1.3? TYPES OF A FUEL CELL .................................................................................... 5?
1.4? ADVANTAGES OF THE PEM FUEL CELL ....................................................... 6?
1.5? BALANCE OF PLANT .......................................................................................... 8?
1.6? APPLICATION TO A QUADRUPED ROBOT .................................................... 9?
2? BACKGROUND .................................................................................................. 10?
3? PEM FUEL CELL STACK .................................................................................. 12?
3.1? PRINCIPLES ........................................................................................................ 14?
3.2? GOVERNING EQUATIONS ............................................................................... 15?
3.2.1? I-V CHARACTERISTIC .......................................................................... 15?
3.2.2? ENERGY BALANCE .............................................................................. 17?
3.2.3? WATER BALANCE ................................................................................ 17?
3.2.4? GDL WITH 2 PHASE EFFECT ............................................................... 20?
3.3? MODEL DIAGRAM ............................................................................................ 23?
4? BALANCE OF PLANT ........................................................................................ 24?
xi?
?
4.1? AIR PROCESSING SYSTEM ............................................................................. 25?
4.1.1? SYSTEM CONFIGURATION AND COMPONENTS ........................... 25?
4.1.2? BLOWER .................................................................................................. 25?
4.2? THERMAL CIRCUIT .......................................................................................... 29?
4.2.1? COOLANT PUMP.................................................................................... 32?
4.2.2? RADIATOR .............................................................................................. 33?
4.2.3? WATER RESERVOIR ............................................................................. 35?
5? ANALYSIS OF A TWO-CELL STACK ............................................................. 36?
6? STACK DESIGN .................................................................................................. 57?
6.1? QUADRUPED ROBOT SPECIFICATION ......................................................... 57?
6.2? STACK SIZING ................................................................................................... 58?
6.3? 20-CELL STACK ANALYSIS FOR CELL EXTENSION ................................. 62?
6.3.1? STACK BEHAVIOR ................................................................................ 63?
6.3.2? DYNAMIC BEHAVIOR .......................................................................... 66?
6.4? 100-CELL STACK ............................................................................................... 71?
6.4.1? MODEL DESCRIPTION ......................................................................... 71?
6.5? DESIGN OF THE AIR PROCESSING SYSTEM ............................................... 74?
6.5.1? SPECIFICATION OF THE AIR PROCESSING SYSTEM .................... 74?
6.5.2? DESIGN OF THE BLOWER ................................................................... 76?
6.5.3? PARAMETRIC STUDY .......................................................................... 81?
6.6? DESIGN OF HEAT EXCHANGER .................................................................... 83?
7? 100-CELL PEM FUEL CELL SYSTEM ............................................................. 90?
7.1.1? STATIC BEHAVIOR ............................................................................... 91?
xii?
?
7.1.2? DYNAMIC BEHAVIOR .......................................................................... 93?
7.1.3? DYNAMIC BEHAVIOR DURING OPERATION MODE ................... 105?
7.1.4? DESIGN PROGRESS ............................................................................. 108?
8? CONCLUSIONS ................................................................................................ 110?
REFERENCES ............................................................................................................... 113?
xiii?
?
LIST OF FIGURES 
Figure 1 Average energy efficiencies of the principal energy conversion systems [1]. ..... 2?
Figure 2 Schematic representation of a PEM fuel cell ....................................................... 4?
Figure 3 A schematic diagram of a PEM fuel cell system [3] ............................................ 8?
Figure 4 Quadruped robot ................................................................................................... 9?
Figure 5 A schematic diagram of a PEM fuel cell stack ................................................... 13?
Figure 6 Mocroscale GDL model ..................................................................................... 21?
Figure 7 Fuel cell stack model schematic block diagram ................................................. 23?
Figures 8 A schematic diagram of balance-of-plant for PEM fuel cell system [3] .......... 24?
Figure 9 A schematic diagram of air processing system .................................................. 25?
Figure 10 A block diagram of blower model .................................................................... 26?
Figure 11 Estimation of blower characteristics and simulation of blower model ............ 29?
Figure 12 A schematic diagram of thermal circuit ........................................................... 30?
Figure 13 Thermal circuit block diagram ......................................................................... 32?
Figure 14 Comparison of the simulated and experimental two-cell stack and individual 
cell voltages ...................................................................................................................... 37?
Figure 15 I-V curve and temperature of a two cell stack under the ramp input of current 
density. .............................................................................................................................. 40?
Figure 16 I-V characteristic, power and overpotential for different operating temperature
........................................................................................................................................... 43?
Figure 17 I-V characteristic, power and overpotential for different stoichiometric ratio . 45?
xiv?
?
Figure 18 I-V characteristic, power and overpotential for different relative humidity .... 47?
Figure 19 I-V characteristic, power and overpotential for different reactant gas pressure49?
Figure 20 I-V characteristic for 9 different cell operating conditions .............................. 50?
Figure 21 Cell power for 9 different cell operating conditions......................................... 50?
Figure 22 Response of a fuel cell to a multi-step current density at different cathode GDL 
thicknesses. (a) Cell voltage, (b) Activation and ohmic overpotential and (c) Membrane 
water content and oxygen concentration near the cathode catalyst layer. ........................ 52?
Figure 23 Dynamic responses of the fuel cells to step changes in current density at 
different cathode GDL porosities. (a) Cell voltage, (b) Activation and ohmic 
overpotential and (c) Membrane water content and oxygen concentration near by the 
cathode catalyst layer. ....................................................................................................... 53?
Figure 24 Dynamic responses of the fuel cells during a step change in current density at 
different membrane thicknesses. (a) Cell voltage, (b) Activation and ohmic overpotential 
and (c) Membrane water content and oxygen concentration near the cathode catalyst layer.
........................................................................................................................................... 55?
Figure 25 I-V characteristic for 4 different geometry parameters .................................... 56?
Figure 26 Cell power for 4 different geometry parameters .............................................. 56?
Figure 27 Stack size and operating voltage
 ...................................................................... 61?
Figure 28 Cell efficiency dependent upon nominal cell potential .................................... 61?
Figure 29 Single cell power, efficiency and current density with selected cell voltage ... 62?
Figure 30 I-V polarization curve of the 20 cell stack ....................................................... 63?
Figure 31 Voltage distribution under the constant current density, 0.6 (A cm
-2
) at 2sec . 64?
Figure 32 Temperature distribution under the constant current density 0.6 (A cm
-2
) ...... 65?
xv?
?
Figure 33 Temperature changes during operation with constant current density ............. 66?
Figure 34 Temperature changes during the start-up ......................................................... 67?
Figure 35 Temperature distribution during the start-up .................................................... 69?
Figure 36 Voltage changes during the start-up ................................................................. 70?
Figure 37 Temperature difference .................................................................................... 71?
Figure 38 A schematic diagram of the 100 cell stack
 ....................................................... 72?
Figure 39 100-cell stack block diagram ............................................................................ 73?
Figure 40 A schematic diagram of blower ........................................................................ 77?
Figure 41 Brake horse power for a different wheel diameter ........................................... 82?
Figure 42 Head loss for a different wheel diameter .......................................................... 82?
Figure 43 Head loss for a same mass flow rate with different hub-tip ratio ..................... 83?
Figure 44 The relationship between coolant temperature difference, flow rate and radiator 
frontal area ........................................................................................................................ 88?
Figure 45 Schematic diagram of PEM fuel cell system .................................................... 90?
Figure 46 I-V polarization characteristic curve of each cell ............................................. 92?
Figure 47 Stack power and voltage ................................................................................... 92?
Figure 48 Temperature changes during the start-up
 ......................................................... 94?
Figure 49 Comparison of temperature distribution during the start-up ............................ 95?
Figure 50 Temperature distribution during the start-up .................................................... 96?
Figure 51 Voltage changes during the start-up ................................................................. 97?
Figure 52 Current density applied..................................................................................... 98?
Figure 53 Transient responses of temperature in the middle cells .................................... 99?
Figure 54 Dynamic response of oxygen concentration to the step input ........................ 100?
xvi?
?
Figure 55 Dynamic response of vapor concentration to the step input ........................... 101?
Figure 56 Transient response of liquid water saturation and water content of membranes 
in the middle cells ........................................................................................................... 102?
Figure 57 Transient responses of activation and ohmic overpotential in the middle cells
......................................................................................................................................... 104?
Figure 58 Transient responses of PEM fuel cell stack voltage ....................................... 104?
Figure 59 Current load under the operation mode .......................................................... 105?
Figure 60 Net power and stack power during the operation mode ................................. 106?
Figure 61 Parasitic power of air supply system and thermal circuit during operation ... 107?
Figure 62 Temperature changes during the operation .................................................... 108?
Figure 63 Schematic diagram and design progress of the PEM fuel cell system ........... 109 
xvii?
?
LIST OF TABLES 
Table 1 Classification of the fuel cell ................................................................................. 6?
Table 2 Three-level orthogonal arrays for 9 different cell operating conditions .............. 49?
Table 3 Two-level orthogonal arrays for 4 different geometry parameters ...................... 56?
Table 4 Operating mode of the quadruped robot .............................................................. 58?
Table 5 Geometric parameters and materials for the design
 ............................................. 72?
Table 6 Design parameter for the blower ......................................................................... 81?
Table 7 Operating condition for the heat exchanger ......................................................... 83?
Table 8 Parameters for heat exchanger ............................................................................. 84?
Table 9 Design progress of the PEM fuel cell system .................................................... 109?
 
  
xviii?
?
NOMENCLATURE 
a Water activity  
A Area (m
2
) 
A
c
 Cross area (m
2
) 
b Parametric coefficient (V) 
c Concentration of species (mol cm
-3
) or safety factor 
C Mass concentration (kg m
-3
) 
CFM Cubic feet per minute (ft
3
min
-1
) 
Cp Specific heat capacity (J kg
-1
 K
-1
) 
d Diameter (m) 
D Diameter (m) 
D
h
 Hydraulic diameter (m) 
D
m
 Diffusion coefficient (m
2
 s
-1
) 
<D
m
> Effective diffusivity (m
2
 s
-1
) 
E Open circuit voltage (V) 
E
electric
 Electrical energy (J) 
g Gravitational acceleration (m s
?2
) 
G Mass flux (kg m
-2
 s
-1
) 
h Convective heat transfer coefficient (W m
-2
K
-1
) 
i Current density (A m
-2
) 
I Current (A) 
j
h 
Colburn factor 
J Rotational inertia (kg m
2
) 
k Conductive heat transfer coefficient (W m
-1
 K
-1
) 
K Gas permeability (m
2
) 
L Length (m) 
L
a
 Air brake power (W) 
L
s
 Shaft horse power 
m Mass (kg)  
M Molecular mass (kg mol
-1
) 
n Number of electrons (2 or 4) 
n
d
 Electro-osmotic drag coefficient 
N Speed (rpm) or molar flux (mol s
-1
 m
-3
) or number of pipe 
N
s 
Specific speed (rpm) 
N
u 
Nusselt number 
p Pressure (atm) 
P Power (W) 
Q Heat (J) 
Q
?
 
Heat transfer rate (W) 
R Universal gas constant (J mol
?1
 K
?1
) 
xix?
?
Re Reynolds number 
R
ele
 Electrical resistance (?) 
R
evap 
Evaporation rate (mol m
-3
 s
-1
) 
s Liquid water saturation ratio 
?S Change in entropy (J mol
-1
K
-1
) 
t Thickness (cm or m) or time (sec) 
T Temperature (K) 
U Overall heat transfer coefficient(W m
-2
 K
-1
) 
v Velocity (m s
-1
) 
V Volume (m
3
) or Voltage (V) 
W Work (J) 
x Molar ratio 
z Number of blades 
 
 
Subscript  
act Activation 
amb Ambient air 
an Anode 
bl Blower 
c Capillary 
ca Cathode 
conc Concentration 
cool Coolant 
d Dynamic 
diff Diffusion effect 
ele Electro osmotic effect 
i Cell number 
in Inlet 
j Layer  
k Reactant or product gases  
l Anode or cathode or liquid 
m Motor 
mass Mass concentration 
mem Membrane 
mid Middle point 
ohm Ohmic 
pl Pole 
rad Radiation 
rct Reactant 
res Reservoir 
rl Relative 
s Stator or static 
sat Saturation 
sou Source 
st Stack 
v Vapor 
xx?
?
  
Greek symbols  
? Porosity 
? Water content 
? Mass density (kg m
-3
) 
null Viscosity (kg m
-1
 s
-1
) 
? Overpotential (V) or efficiency 
? Constant parametric coefficient (V) 
null Stress (Pa) 
?
mem
 Membrane conductivity (?
-1
 cm
-1
) 
null Specific mass of gas or volumetric condensation rate (s
-1
) 
null Contact angle (
o
) 
? Flux linkage (V s rad
-1
) 
null Angular velocity (rad s
-1
) 
 
 
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
 
1?
?
1 INTRODUCTION 
?
1.1 THE FUEL CELL 
Currently, research about renewable alternative energy is energetically pursued to 
solve economical issues from highly priced gas to environmental issues like global 
warming. In other words, research about renewable alternative energy is not only for 
economic benefit, but also for the preservation of human beings? very existence. 
Conventional fossil fuels which meet 80% of the energy demands have two major 
disadvantages. The first problem is that the amount of fossil fuels is limited. Secondly, 
fossil fuels cause severe pollution which lead to global warming, weather changes, the 
Antarctic icecap melting, an increasing sea level, acid rain and destruction of soil and 
forests. Therefore, recently, many industries and governments have concentrated on 
researching and developing renewable alternative energy. 
2?
?
?
Figure 1 Average energy efficiencies of the principal energy conversion systems [1]. 
 
There are many technologies being developed to replace conventional energy 
sources. Specifically, hydrogen is one of the most promising energy fuels of the future. A 
hydrogen energy system was actually proposed in early 1970?s. However, its research has 
been drastically increased over the last 10 years by many universities and research 
institutes. Hydrogen has a lot of advantages over other contenders. It is the most 
lightweight fuel as well as the most efficient and clean energy source [2]. 
A fuel cell is an electrochemical energy conversion device which converts the 
chemical energy of fuel such as hydrogen into electrical energy like DC electricity. One 
of the advantages of a fuel cell is that it can efficiently convert hydrogen into electricity 
by electrochemical processes. This conversion efficiency is much higher than other 
traditional energy conversion devices such as steam engines to turn generators. Several 
3?
?
steps were required to make electrical energy in the past. For example, heat, produced 
through burning fuel, is used for making a working fluid like water vapor. The vapor is 
used to rotate a turbine to make mechanical energy which is converted to electrical 
energy by use of an electric generator. However, a fuel cell is able to minimize this 
process of making electricity without the use of any mechanical device, and this is the 
reason why the efficiency of fuel cell is usually higher than that of traditional internal 
combustion engines. Moreover, it guarantees a fuel cell a longer life than other 
mechanical energy conversion devices in terms of mechanical endurance because it has 
no moving parts like pistons, shafts, or bearings. 
Fuel cells also have big advantages in environmental issues. An impurity like 
sulfur causes performance and endurance of fuel cells to decrease largely. Therefore, 
impurities are clearly removed before it is delivered and used in a fuel. Consequently, it 
makes clean exhaust gas without any SO
x
. Moreover, pure water from the chemical 
reaction between hydrogen and oxygen is a pivotal factor in making the fuel cell 
environmentally friendly. Thus, fuel cells are highly efficient, environmentally friendly 
energy conversion devices that can potentially replace conventional internal combustion 
engine. 
 
1.2 PRINCIPLES OF A FUEL CELL OPERATION 
A fuel cell has an oxidizing agent on the catalyst layer where a chemical reaction 
occurs with fuel and oxygen. Hydrogen is usually used as a fuel because its 
electrochemical reaction speed is the fastest among the fuel cell fuels. Oxygen is used as 
an oxidant, which abundantly exist in the air and can be used to supply oxygen. When 
4?
?
hydrogen is used as a fuel and oxygen is used as an oxidant, fuel cells are one of the 
cleanest energy conversion devices, only producing pure water and heat as the 
byproducts.  
H
+
?+
+? eHH
Anode
22
2
EnergyHeat
EnergyElectrical
lOHOH
+
+
?+ )(
2
1
ReactionOverall
222
O
Membrane Electrode Assembly
A
node C
a
t
a
ly
st
C
a
t
hode C
a
taly
st
Polymer 
Electrolyte 
Membrane
G
a
s D
i
f
f
u
s
io
n 
Lay
er
G
a
s D
i
f
f
u
s
io
n 
Lay
er
H
+
H
+
H
+
H
H
O
O
H
+
H
+
-
-
-
-
-
-
-
H
H
H
H
H
+
H
+
22
442
Cathode
OH e HO
+
++?
H
H
O
v
-
v
+
?
Figure 2 Schematic representation of a PEM fuel cell 
?
 
A fuel cell causes an electrochemical reaction to occur between both faces of a 
permeable membrane which transmit no gas but only specific ion selectively. On one side, 
reaction occurs, which produces specific ions, but on the other side, a reaction occurs, 
which consumes specific ions.  
For example, as shown in the figure above, electrons generating reaction occurs as 
the hydrogen being supplied is converted into hydrogen ions. The hydrogen ions 
5?
?
produced are transported through a membrane to the cathode catalysts in the fuel cell, 
while the electrons produced move to the opposite side through the outer transport circuit, 
reacting with oxygen supplied and then produces water. In sum, each oxidation reaction 
of hydrogen and reduction reaction of oxygen occurs on both sides of the membrane as 
water and electricity are produced.  
Fuel cells use platinum as a catalyst that exists on both sides of the membrane and 
allows for electrochemical reactions. A catalyst layer that separates electrons from the 
oxidation reaction of the fuel is called the anode. The other catalyst layer which 
consumes electrons from the reduction reaction with the oxygen supplied is called the 
cathode. Therefore, about 50% of energy in fuel is converted into electricity, not by 
combustion, but by an electrochemical reaction between the fuel and air. The remainder 
is converted to heat. 
 
1.3 TYPES OF A FUEL CELL 
Fuel cells can be categorized by the operating temperature and fuel delivery 
methods, but can be more specifically classified by the type of electrolyte used. Polymer 
electrolyte membrane fuel cells (PEMFC) have some specific advantages over other fuel 
cells. They can reduce thermal loss by operating at lower temperatures, and their smaller 
size is more appropriate for automobiles and portable devices. This is the reason why 
most automotive companies research polymer electrolyte membrane fuel cells and why 
almost all the fuel cells used in automobiles are polymer electrolyte membrane fuel cells.  
 
 
6?
?
 
Classification criterion 
Operating temperature 
(?C) 
Electrolyte 
Power 
density 
(kW m
-2
) 
Electrical 
efficiency 
of stack 
(%) 
Low 
temperature 
60-100 AFC (Alkaline fuel cell) 0.7-8.1 45-60 
50-100 
PEMFC (Polymer 
electrolyte membrane fuel 
cell) 
3.8-6.5 40-55 
50-200 
DMFC (Direct methanol 
fuel cell) 
- 40 
160-210 
PAFC (Phosphoric acid fuel 
cell) 
0.8-1.9 40-50 
High 
temperature 
650-800 
MCFC (Molten carbonate 
fuel cell) 
0.1-1.5 
50-60 
(system) 
800-
1000 
SOFC (Solid oxide fuel 
cell) 
1.5-2.6 50-65 
Table 1 Classification of the fuel cell 
?
1.4 ADVANTAGES OF THE PEM FUEL CELL 
As mentioned before, fuel cells have a higher efficiency than traditional internal 
combustion engines, so they are very suitable for automobile propulsion system. For 
instance, the electrical efficiency of polymer electrolyte membrane fuel cells is higher 
7?
?
than that of a comparably sized internal combustion engine (about 55% optimistically 
compared to 37.6% at full load [1]. Moreover, fuel cells which use hydrogen as a fuel 
exhaust harmless gases like air and water vapor. They are also very resistant to corrosion 
because the reaction product is not a liquid electrolyte, but only water.  
Polymer electrolyte membrane fuel cells operates at a low temperature, that is, 
about 80?C, which make polymer electrolyte fuel cell systems stable and suitable for use 
in automobiles and portable devices because of a smaller possibility of accidents than 
conventional internal combustion engines. The low operating temperature allows polymer 
electrolyte membrane fuel cells to reduce start-up time much more than other types of 
fuel cells. These are the main advantages of the polymer electrolyte membrane fuel cell 
over traditional energy conversion devices. Furthermore, it can be modularized so that the 
amount of power generation is increased easily by attaching additional fuel cell stacks to 
meet demands. This feature enables us not only to gradually increase the size of energy 
conversion devices but also to minimize initial cost for energy conversion devices. 
Having no moving parts also allows a fuel cell to increase its life span. Currently, the life 
span of a fuel cell is about 4,000 hours which should be further increased to 40,000 for 
use in power plant stations. Its silent operation is another advantage for application in 
portable devices and military devices. However, polymer electrolyte membrane fuel cells 
should be hydrated properly with a supply of humidified air to maintain appropriate 
hydrogen ion conductivity. Despite some disadvantages, polymer electrolyte membrane 
fuel cells are the best potential energy conversion device for automobiles. 
 
8?
?
1.5 BALANCE OF PLANT 
Auxiliary components such as the air supply system, fuel delivery system, and 
thermal system including a coolant pump, radiator, and fan are essential for the 
successful operation of a fuel cell system. The air supply system and fuel delivery system 
supply reactant gas, oxygen and hydrogen to the stack. A humidifier in the air supply 
system is another important part because a well hydrated membrane can transport 
positive ions with less ohmic loss. Therefore, the reactants of both the anode and cathode 
should be humidified appropriately. A thermal system manages heat produced from the 
cell stack and maintains the desired operating temperature which is pivotal for the 
performance and durability of a fuel cell.  
Hydrogen 
Tank
M
B
Air
Blower
PEMFC
Radiator & Fan
Inlet Manifold
Outlet Manifold
Inlet Manifold
Humidifier
Outlet Manifold
Bypass Valve
S
S
S
Pump
Water 
Reservoir
?
Figure 3 A schematic diagram of a PEM fuel cell system [3] 
?
9?
?
These auxiliary components constitute the balance of plant (B.O.P). They work 
dependently of each other because they are connected and influence each other 
thermodynamically and materially in order to accomplish optimum thermal efficiency. 
Therefore, each system is made to be able to engage with other components during 
operation. This balance of plant costs about 70% of the entire fuel cell system. 
Consequently, the development of the balance of plant, as well as the fuel cell stacks, 
which needs to satisfy technical specification, safety and economical demand are most 
important. 
 
1.6 APPLICATION TO A QUADRUPED ROBOT 
?
Figure 4 Quadruped robot 
?
One of the main motivations for developing fuel cell based propulsion system for a 
multi-legged type robot is the possibility of the quietness and compact size of the purely 
electric propulsion systems, namely high energy and power density. Fuel cell propulsion 
systems can also profit from various possibilities for improving the fuel economy with 
respect to ICE-based propulsion systems [2]. 
10?
?
2 BACKGROUND 
PEM fuel cell is an electrochemical energy conversion device whose durability 
and performance is easily affected by operating condition [2]. Moreover, heat and water 
is byproduct during an operation of the PEM fuel cell and its appropriate management 
plays an important role in maintaining and improving the performance of the fuel cell. 
Therefore, PEM fuel cell system should be operated with careful water and heat 
management for the performance and durability.  
In addition, PEM fuel cell is stacked with many components in the cell such as 
end plate, gas channel, coolant channel, gas diffusion layer and membrane. Each physical 
properties like pressure, humidity and temperature is hardly measured for the analysis 
because it is sealed fully in the enclosure of the fuel cell. Therefore, modeling fuel cell 
can offer better information and knowledge about physical characteristic of the fuel cell 
stack. 
There are ways to model PEM fuel cell. One is modeling using high order 
equation for two dimensional or three dimensional model using CFD techniques. C.Y. 
Wang proposed PEM fuel cell model to describe water generation and liquid water phase 
in the GDL [4]. It provides information about physical phenomena which occur in the 
stack. However, it has a limitation to explain behavior in the multi-cell stack due to the 
slow computational time and over abundance of the parameter and variables used in 
computer simulation. Other way is using empirical model as low order. J.C. Amphlett 
proposed lumped thermal model which can describe heat generation and
11?
?
 voltage characteristics [5]. The model can be simple but ignore many physical 
phenomena such as the coupling of fluid and thermo dynamics. Therefore, there is need 
for combination of each advantage and minimization of disadvantage. Shan [6] has 
proposed model which consider effects of a varying load n dynamics of temperature. It is 
more predictive than empirical model and less complicated than analytical model using 
medium order equation. However, the model which is mentioned above is developed with 
single cell or small number of cells in which the temperature effect can be neglected. In 
the real application such as mobile devices and automobiles, the PEM fuel cell stack is 
developed as a multi cell stack by connecting a lot of cells about 400 in series. Therefore, 
the model with large number of cell is needed to be developed to explain the interaction 
of temperature effect and dynamic behavior of physical properties at varying loads. 
Based on analysis including comparisons a various operating condition and geometry 
factor in the model, the model?s operating condition is optimized and the 100-cell stack is 
developed based on 20-cell stack analysis. The operation of a 100-cell stack under 
varying conditions, and the transient behavior of a step load are analyzed. Temperature 
distribution and the change of overpotentials along with the associated cell voltages were 
the major parameters investigated.  
12?
?
3 PEM FUEL CELL STACK 
In a fuel cell stack, there are complex correlated behaviors between 
thermodynamics and fluid mechanics. These complex phenomena have been explored 
either by empirical methods using curve fitted equations or numerical methods exploiting 
advanced two- or three- dimensional computational fluid dynamics (CFD) techniques. 
The empirical method is good due to less calculating time and less memory consumed, 
but is usually valid for a specific steady state. Therefore, empirical methods were not 
used in analysis of dynamic behavior.  
While the CFD techniques can explain dynamic behavior in more detail, the 
techniques have been limited to analyzing a single cell because of the enormous 
computational time required for a higher number of cells due to a lot of grid points and 
the corresponding parameters that need to be determined for simulations. Consequently, a 
dynamic model for the PEM fuel cell stack is proposed considering temperature 
distribution and a two-phase effect using empirical equations. The model was used to 
analyze the dynamic behavior of cell stack.  
A fuel cell is an electrochemical energy conversion device. A PEM fuel cell is 
composed of various layers including bipolar plates, gas diffusion layers, membrane 
electrode assembly (MEA), and catalyst layers. The model is constructed using two cells 
connected to each other.  
 
13?
?
 
?
1 Anode end plate 6 
Anode catalyst, MEA, 
Cathode catalyst  
2 
Anode bus plate (Current 
collector) 
7 Cathode gas channel plate  
3 Anode cooling plate 8 Cathode cooling plate 
4 Anode gas channel plate 9 
Cathode bus plate (Current 
collector) 
5 Gas diffusion layer 10 Cathode end plate 
Figure 5 A schematic diagram of a PEM fuel cell stack 
?
  
Some assumptions for the modeling are made. All gases in the model obey the 
ideal gas law. There is no pressure gradient between the anode and cathode side which 
means no convection, only diffusion when gas transport is considered. There is no gas 
pressure drop from the inlet to the outlet of the gas channel. The temperature gradient is 
linear across the layers in a fuel cell. The thermal conductivity for the materials in a fuel 
cell is constant. Anodic overpotential is negligible. There is no current density gradient 
across the cathode catalyst layer which means that the reactants completely reacted as 
14?
?
soon as they reach the cathode catalyst layer surface. And finally, latent heat during a 
phase change is not considered. 
Based on these assumptions, the model is developed which is described in the 
following sections. Static behaviors of a cell model are based on voltage equations and 
the effects of temperature distribution in a cell, water balance in the membrane, and gas 
dynamics in the cathode gas diffusion layer with a two-phase phenomenon. 
 
3.1 PRINCIPLES 
In the system, there are electrochemical reactions at both the anode and cathode 
side. 
Anode: 
2
22HHe
+?
?+
 
???????
(3- 1) 
Cathode: 
22
1
22
2
OH eHO
+?
++?
? ?????
(3- 2) 
To be specific, these reactions occur at the interface between a conductive 
electrolyte and a conductive electrode. The catalyst usually has a porous structure 
because the gas reactants should be able to reach the reaction location. The product from 
this reaction is water that should be ejected properly. The reaction between hydrogen and 
oxygen usually produces not only water, but also electricity and reaction heat. 
222
1
2
electric
HOHOQE+?++
???????
(3- 3) 
?
      
 
 
15?
?
3.2 GOVERNING EQUATIONS 
The stack model is developed as a single cell at first using some empirical 
equations. A single cell consists of layers such as the bipolar plate, membrane, catalyst 
and gas diffusion layer. 
 
3.2.1 I-V CHARACTERISTIC 
A voltage is calculated by subtracting activation overpotential, ohmic 
overpotential through membrane and concentration overpotential at cathode catalyst layer 
from open circuit voltage E
i
. 
coniohmiactiii
EV
,,,
??? ???=
???????
(3- 4) 
Open circuit voltage is calculated by considering the relationship between 
chemical and electrical potential energy [7]. 
() ()
?
?
?
?
?
?
+??+???=
??
22
,,
53
ln
2
1
ln103085.4)15.298(1085.0229.1
OiHiiii
ppTTE
?
(3- 5) 
As shown in the Butler-Volmer equation, there should be a difference between 
electrode potential and equilibrium potential in order to accomplish the electrochemical 
reaction and produce current. This is activation polarization which is related to the slow 
reaction speed. Therefore, activation polarization loss becomes less as current density 
increases. The reduction reaction of oxygen usually needs more activation overpotential 
than the oxidation reaction of hydrogen, which means that reduction reaction is much 
slower than oxidation reaction [7]. 
)ln()ln(
4,321,
2
ITcTT
iOiiiacti
??+??+?+= ?????
?? ? ?
(3- 6) 
16?
?
2
2
,
,
6
5.08 10 exp( 498 / )
iO
iO
p
c
T
=
??
???
????
(3- 7) 
where ?
i
 
is constant parametric coefficients, I is current (A), 
2
O
c is oxygen concentration. 
A resistance loss occurs due to the resistance caused by ion transport at interfaces 
such as electrolyte membranes and electron transport near the electro-conductive 
components like electrode surfaces in fuel cells. These losses can be explained by 
considering Ohm?s law [8]. 
memi
mem
ohmiohmi
ti
Ri
,
,,
?
?
?
=?=
?? ?
?
????
(3- 8) 
where ?
ohm
 is ohmic overpotential, R
ohm
 is the ohmic resistance, t
mem
 and ?
mem
 are the 
membrane thickness and conductivity. The membrane conductivity is dependent upon 
temperature and water content [8][9]. 
,
,
11
exp 1268
303
0.005139 0.00326
imem i
i
iime
T
??
??
??
??
=? ?
????
??
??
=??
??
?
????
(3- 9) 
where ?
mem
 is membrane water content. 
There is a concentration gradient near an electrode due to the electrochemical 
reaction as mentioned earlier. The concentration gradient is the cause of concentration 
overpotential [10]. 
max
,
max max
ln ln 1
1
icon
iRT i
b
nF i i
?
?? ??
==???
?? ??
?
?? ??
?
?
????
(3- 10) 
where ?
con
 is concentration overpotential, b is parametric coefficient (V), I is actual 
current density (A cm
-2
), I
max 
is maximum current density (A cm
-2
). The voltage of the 
stack is the sum of the voltage of each single cell in the stack. 
17?
?
stack i
i
VV=
?
?? ? ?
?
?
????
(3- 11) 
Where, V
stack
 
is the voltage of the stack and V
cell
 
is the voltage of a single cell. 
 
3.2.2 ENERGY BALANCE 
Energy flow inside of the stack consists of thermal conduction between layers, 
mass flow into each volume, and convection from reactant flow and coolant under the 
assumptions that physical properties are independent of temperature and are isotropic. In 
other words, the energy exchange at the interface of the control volume and the internal 
heat source constitutes the energy balance. This energy flow can be explained by 
examining first law of thermodynamics [11]. 
()
j
j j k k in j convection conduction sou
j
mass flow
dT
mCp mCpT T Q Q Q
dt
??= ? ?+ + +
??


	

?? ?
(3- 12) 
Where, the first term presents mass flow, the second and third terms are convection, the 
fourth term is conduction heat transfer, and last term is the heat source term. 
Entropy change, activation overpotential and ohmic overpotential cause heat 
generation during the reaction [12]. 
i
sou i ele
TS
QiA iR
nF
?
????
=? ? + +?
??
?
??

?? ? ?
???
(3- 13) 
?
3.2.3 WATER BALANCE 
A membrane should offer high conductivity appropriately and be chemically and 
mechanically stable under the fuel cell?s operating conditions. The water balance in the 
18?
?
membrane is managed by the electro-osmotic driving force caused by the difference in 
electrochemical potential at both the anode and cathode, while the diffusion caused by the 
difference in concentration results in a concentration gradient. 
Conductivity of a polymer electrolyte membrane is mainly dependent on the type 
of structure and water content in a membrane. The maximum water content of a 
membrane can be determined by the phase of water in the membrane. Usually, the 
forming of liquid water increases water content in the membrane.  
The water content in the membrane is calculated with consideration of mass flow 
at the boundaries within the membrane [11]. 
()
22
22
2
,,
,
,
,,
,,
, ,,, ,,, ,,, ,,,
/
0.0126 /
iHOmass HO
imem
dry mem
iHOmass HO
mem
iHOmass mem
i mem i ele mem an i ele mem ca i diff mem an i diff mem ca
CM
CM
M
dC At
mm
dt
?
?
=
??
??
==?++
?
(3- 14) 
where C is mass concentration (kg m
-3
), M is mole mass (kg mol
-1
), ? is membrane dry 
density and A
cell
 is fuel cell area (m
2
). 
As mentioned earlier, there are different ways to transport the water. The electro-
osmotic driving force is one of the causes that induce water mass flow across the 
membrane layer. The water is produced from the electrochemical reaction while the rate 
of water production is a function of load [8]. 
2
,
2
H O generation
i
N
F
=
?????????
(3- 15) 
where N
H2O,generation
 
is the rate of water generation (mol cm
-2
 s
-1
). 
19?
?
Water is carried by the electro-osmotic drag force with the protons which move 
from anode to cathode. The electro-osmotic drag coefficient is a function of water content 
[8]. 
,, , ,,ielememl water dil
i
mMAn
F
=??? ?
??????
(3- 16) 
19
,
2
,,,
104.305.00029.0
?
???+?=
lililid
n ??
?? ? ? ? ?
(3- 17) 
 The concentration gradient at the boundaries of the membrane layer increases due 
to water generation and water transport caused by the electro-osmotic drag force. 
Diffusion, W
diff,mem
, is caused by this concentration gradient. In some cases, the 
concentration gradient is capable of reverse diffusion, from cathode to anode, during 
operation [8]. 
( )
,,
,, , ,
il imid
i diff mem l water water i l
mem
CC
mMAD
t
?
=??
?? ? ? ? ?
(3- 18) 
memi
mem
memdry
midili
mem
memdry
li
M
C
M
C
,
,
,,
,
,
, ?
?
?
?
?=?=
?
?????
(3- 19) 
The diffusion coefficient is proposed in the form of an empirical equation [8]. 
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
??=
imem
liliwater
T
DD
,
,,,,
1
303
1
2416exp
?
?
?????
(3- 20) 
()()
()
5.4
35.4
23
2
1025.1
367.1310
32110
10
,
,
,
,
6
.
6
.
6
6
,,
?
>>
??
>
?
?
?
?
?
?
?
?
??
?+?
=
?
?
?
?
li
li
li
li
li
li
li
D
?
?
?
?
?
?
?
?? ? ? ?
(3- 21) 
The water content is also calculated as a function of water activity [8]. 
20?
?
()
3
13
01
8.16
14.114
3685.3981.17043.0
,
,
,
,
3
,
2
,,
,
?
>?
>?
?
?
?
?
?
?+
?+???+
=
li
li
li
li
lilili
li
a
a
a
a
aaa
?
?? ? ?
(3- 22) 
The water activity is a relationship between water vapor partial pressure and 
saturation pressure [8]. 
lisat
liv
li
P
P
a
,,
,,
,
=
??????
????
(3- 23) 
?
3.2.4 GDL WITH 2 PHASE EFFECT 
Water in the form of a liquid plays a pivotal role in improving performance by 
increasing proton conductivity or decreasing performance by impeding the transport of 
gas reactants. Not only the gas phase, but also the liquid phase - which has not been 
previously considered - are taken into account, diffusivity should be modified to effective 
diffusivity which represents the diffusion characteristic of liquid water and water vapor in 
a capillary in order to show the dynamic behavior of the water. It is a function of the 
diffusion coefficient, porosity of diffusion layer, and liquid water saturation ratio [13]. 
()
vaporOm
V
kV
ks
V
V
ksDkD
pore
liquid
GDL
pore
mm
,
)(
)(,
)(1
11.01
11.0
)(
2
2
785.0
=
==
?
?
?
?
?
?
?
?
?
=??
?
?
?
??
????
(3- 24) 
where <D
j
> is effective diffusivity (m
2 
s
-1
), D
j
 is diffusion coefficient (m
2
 s
-1
) at a single 
phase, ? is porosity of diffusion layer (-), s is liquid water saturation ratio (-), V
p
 is pore 
volume of the GDL (m
3
), V is total volume of GDL(m
3
), and V
l
 is volume of liquid water 
(m
3
). 
21?
?
 Effective diffusivity of gas diffusion layer is developed by using fiber screen 
model [13]. Real fiber structure is shown as Figure 6 (c). Its randomly overlapped 
structure is hard to be modeled with pore space. Moreover, the GDL used in PEM fuel 
cell is stacked fibrous structure. Therefore, it is hard to define and measure pore portion 
inside the GDL. The microscale model for GDL is developed based on the assumption 
that the pore spaces is made by four interwoven fiber as a square shape. This makes 
modeling GDL much easier and its porosity can be calculated as a volumetric ratio. 
However, its height is assumed to be a same size. Consequently, porosity of the GDL can 
be derived as below [13]. 
??
(a)?A?plane?of?screen?made?of?fibers??? (b)?Stacked?fiber?screen?model?
?? ? ?
(c) ?Random?fiber?structure?? ? (d)?Unit?square?pore?space?
Figure 6 Mocroscale GDL model 
?
22?
?
2
2
()
p
pf
d
dd
? =
+
   
??
????
(3- 25) 
where d
p 
is pore size and d
f
 is the representative fiber diameter (m). 
 The mass conservation is in the control volume, that is, each section in the gas 
diffusion layer is considered in order to describe the volume of liquid water. There are 
two types of behavior in the control volume, liquid water flow by concentration gradient 
and evaporation [13]. 
,,
evap p v
l
llinlout
RVM
dV
mm
dt
?
?
=? ?
 
??
????
(3- 26) 
where ?
l
 is density of liquid water (kg m
-3
), 
l
m is mass flow rate of liquid water (kg s
-1
), 
R
evap
 is molar evaporation rate (mol sec
-1
 m
-3
), M
v
 is molecular weight of vapor (kg mol
-1
) 
and ? is porosity of the diffusion layer (-).  
 Evaporation rate, which is related to the volume of liquid water, is mainly 
dependent on the difference between the saturation pressure and water vapor pressure 
[13]. 
RTcp
RT
kpTp
kR
vv
vsatv
evap
=
?
= ,
)()(
)(
,
?
??
????
(3- 27) 
where ? is volumetric condensation coefficient (s
-1
), P
v,sat 
is saturation pressure (Pa), P
v
 is 
partial pressure of vapor (Pa) in each region, c
v
 is molar concentration of vapor (mol m
-3
), 
R is ideal gas constant (J kg
-1
 K
-1
) and T is absolute temperature (K). 
?
23?
?
3.3 MODEL DIAGRAM 
The structure of the model is developed by the SIMULINK, which has a graphical 
interface and block diagram. The program allows us to analyze dynamic behavior 
arranging theses icons into a block diagram representation of a dynamic process.  
As explained above, each layer is connected with consideration of heat and mass 
transfer. The layers are developed using continuity equation and heat transfer equation at 
the end plate, coolant channel and gas channel. Specifically, catalyst layer uses 
electrochemical equation such as open circuit voltage, activation overpotential and 
concentration overpotential equation. Membrane uses important equation which decides 
proton conductivity, membrane water content and ohmic overpotential. Finally, GDL of 
cathode is developed using a two-phase effect for the explanation of flooding phenomena. 
Like this, anode and cathode part are developed and composed. 
inca
m
,

i
inan
m
,

Gas diffusion layer
Cell voltage
GDLgas
Q
?

Gas diffusion layer
GDLgas
Q
?

catlGDL
Q
?

gascoolant
Q
?

gascoolant
Q
?

Gas flow channel
catlGDL
Q
?

Membrane
memGDLwater
m
?,

memGDLwater
m
?,

mem
?
GDLgaswater
m
?,

GDLgaswater
m
?,

outca
m
,

outan
m
,

generatedwater
m
,

catlGDLO
m
?,
2

Cathode
Anode
cell
V
Gas flow channel
GDLgasO
m
?,
2

GDLgasH
m
?,
2

Catalyst layer
Catalyst layer
catlGDLH
m
?,
2

+
H
+
H
2
O
p
2
H
p
catl
T
?
Figure 7 Fuel cell stack model schematic block diagram 
24?
?
4 BALANCE OF PLANT 
The fuel cell system is composed of many subsystems such as an air processing 
system, a fuel processing system, thermal management and water management as well as 
the fuel cell stack system. There are also many components to these subsystems such as a 
blower, hydrogen tank, humidifier, water reservoir, purge valve, and a deionized water 
pack. Proper operation of each subsystem and its components makes a successful fuel cell 
system work. Fuel cell system is organized like the following. 
Hydrogen 
Tank
S
M
BL
Air
Blower
PEMFC
Radiator & Fan
Inlet Manifold
Outlet Manifold
Inlet Manifold
Water 
Reservoir
Humidifier
Pump
Pump
Condenser
Condenser
Outlet Manifold
Thermostat
?
Figures 8 A schematic diagram of balance-of-plant for PEM fuel cell system [3]? 
?
25?
?
4.1 AIR PROCESSING SYSTEM 
?
4.1.1 SYSTEM CONFIGURATION AND COMPONENTS 
?Air is one of the important reactants needed in order to produce electricity in a 
fuel cell. Considering the stoichiometric ratio, the proper amounts of fuel and air should 
be supplied. For this, the air supply system is composed of a blower driven by a motor 
humidifier, manifold and others.  
The blower serves to provide the specific amount of air required. Moreover, a 
specially designed blower, such as a high pressure blower, is preferred for a fuel cell 
system. A brushless motor is usually used in fuel cell systems to prevent leaking fuel 
from being ignited and exploded. 
?
Figure 9 A schematic diagram of air processing system 
 
4.1.2  BLOWER 
?
The current blower model was originated from a characteristic curve which was 
composed of the efficiency and flow rate curves depending on head loss parameters. 
26?
?
?
dt
J
bl
1
bl
V
+
?
mbl
plmbl
N
,
,.
2
?? ?
?
?
?
?
?
?
?
? ?
?
?
?
?
?
?
?
?
?
?
mbls
mblplmblmbl
R
N
,,
,,,,
22
3 ?
?
m
?
?
+
bl
?
bl
?
ca
p
?
+
amb
p
?
?
( )uf
( )ug
bl
m
bl
p?
?
ambbl
mbl
??
?
?
,
bl
?
*
p
bl
?
?
Figure 10 A block diagram of blower model 
?
The blower is operated by a DC brushless motor, but its electrical characteristics 
are not considered. The output of the blower model is significantly affected by the inertia 
of the rotor, the required torque and the applied torque. 
 
4.1.2.1  MOTOR-INPUT VOLTAGE 
?
27?
?
4.1.2.1.1  BLOWER SPEED 
The dynamic behavior of the blower speed is described by the moment of the 
inertia. The blower torque is the difference between the motor toque applied and the 
blower?s load torque. 
blmblbl
J ??? ?=
         (4- 1) 
Where ?
m
 is motor torque (N m) and ?
bl
 is blower torque (N m). 
The motor torque is a function of the power input and motor speed. The required 
blower torque is also a function of power and speed but some additional losses occur due 
to the blower?s loading torque. Thus the equations for torque are 
?
?
?
?
bl
bl
m
m
P
P
=
=
          (4- 2) 
where the required power of the blower, P
bl
 (W), is calculated by the thermodynamic 
power equation. 
The motor torque is calculated and linearized under the three phase balanced 
operation representation [5]. 
28?
?
,
,
1
bl bl m
bl
bl m
sbl
W
d
dt J
?
?
?
?
??
=?
??
??

        (4- 3) 
bl bl
bl
bl amb
mp
W
??
?
=


          (4- ) 
,, ,,,
,, , ,
,,
3
22 2
bl m pl bl m plbl m
bls m bl m rqs bl m bl m bl
sblm
NN
V
R
? ??
??? ??
?? ??
=??
??? ?
?? ??
?? ??
? ???
   (4- 5) 
where ?
bl,m
 is a torque produced by the motor (J), W
bl
 is blower work (J), R
s,bl,m
 is a 
function of the stator resistance(Ohm), ?
bl,m
 is flux linkage (V s 
-1
 rad
-1
) and n
bl,m,pl
 is the 
number of the poles with the stator voltage, V
bl,m 
(V). 
?
4.1.2.2  AIR FLOW RATE 
The total mass flow rate is calculated as a function of the angular velocity and 
pressure, and the efficiency given as a function of the flow rate and the angular velocity 
[16]. The blower map is used to determine the air flow rate through the blower [14]. 
()
( )
()
2
*3*5*422
*3*422
10.29 0.72 10 2.07 10 , 9 10 /
0.9 0.8210 , 910 /
bl bl
bl bl
mpppPasrad
m p p Pa s rad
?
?
???
??
=?? ? ? ? ?+ ? ?? ?
=???+ ? ?? ?


  (4- 6) 
Where p* is 
?
?
?
?
?
?
?
? ?
2
bl
ambca
pp
?
.
 
The parameters of the blower were extracted from desired specifications and 
characteristic data delivered by PADT (Phoenix Analysis & Design Technologies) [15]. 
In a low pressure system, 
29?
?
1
manifold
atm
p
p
??
?
??
??
, 
manifold bl atm
TTT??        (4- 7) 
The air flow characteristics are air mass flow, temperature, pressure, relative 
humidity, and the oxygen mole fraction. But in a low pressure system, temperature, 
pressure, and relative humidity are the same as atmospheric ones while the oxygen mole 
fraction is 0.21. 
0 2000 4000 6000 8000 10000 12000 14000 16000 18000 20000
0.0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1.0
1.1
1.2
 Model
 Estimation
 
 
Flo
w
 R
a
te
 (C
M
M
)
RPM
?
Figure 11 Estimation of blower characteristics and simulation of blower model 
?
4.2 THERMAL CIRCUIT 
Products from the chemical reaction in a PEM fuel cell include not only 
electricity and water but also heat. Excessive heat is able to raise the temperature of a 
PEM fuel cell abnormally causing damage to the stack. An appropriate temperature is 
30?
?
critical in order to increase the performance of a PEM fuel cell, such as conductivity in 
the membrane which is directly related to the water content in the membrane. Therefore, 
an appropriate thermal management system plays an important role in the successful 
operation of fuel cell stacks.  
?
Figure 12 A schematic diagram of thermal circuit 
?
The thermal circuit is composed of a radiator, fan, coolant pump and water 
reservoir. The selection of the proper size of a pump and a radiator is necessary because 
the thermal system should be able to reject the excessive heat generated in the PEM fuel 
cell. An approximation of the heat produced can be calculated as a function of current, 
entropy, overpotential and ohmic overpotential [12]. 
31?
?
i
sou i ele
TS
QiA iR
nF
?
????
=? ? + +?
??
?
??

     
  (4- 8) 
 The cell stack generates approximately 20 kW of heat during maximum output.
?
Heat produced can be assumed to be the activation and ohmic voltage loss which varies 
depending on load current density. 
()
ohm act conc
VE V V V V=?? ?? ??
     
  (4- 9) 
Therefore, the capacity of a coolant pump and a radiator should be decided based 
on the operating conditions and required specifications. The maximum possible heat 
produced will occur at maximum current density within range of operation, and the 
thermal components should satisfy all of the cooling demands.  
Heat rejected by the thermal system from the stack to coolant is determined by the 
coolant flow rate and the intended temperature difference between the stack inlet and 
outlet.  
,
()
coolant generated transfer coolant p coolant coolant
QQ mC Tm?=?=???


  
  (4- 10) 
After the heat produced and the required temperature difference is calculated, the 
necessary coolant flow rate can be calculated. Heat is needed to reject about 20 kW from 
the calculation using equation above for heat generation rate during operation of the PEM 
fuel cell at maximum power.  
32?
?
coolant
m
coolant
m
amb
T
incoolant
T
,
?
+
air
m
()uf
radiator
h
coolant
frontrad
Cp
A
_
5.0 ?
? ?
+
outcoolant
T
,
( )ug
fan
P
res
T
Ah ?
amb
T
?
+
?
+
c
Cp
?
?
+
res
T
outcoolant
T
,
resres
Cpm
t
?
?
?
Radiator & fan
Water reservoir
stack
T
?
+
Control
T?
coolant
T
Coolant pump
 
Figure 13 Thermal circuit block diagram 
 
4.2.1 COOLANT PUMP 
 The coolant pump plays a pivotal role in circulating the coolant in the circuit from 
the heat source to the heat exchanger. Under the assumption that heat produced during 
the operation is rejected fully by the coolant, the required mass flow rate is calculated 
using the calculated heat generation, type of coolant, i.e. heat capacity of coolant, and 
temperature difference between inlet and outlet. 
33?
?
_stack gen
cool
cool
Q
m
Cp T
=
??


     
  
  (4- 11) 
 There are some different ways to calculate the heat balance. Generally, the 
chemical energy of the fuel is converted into electricity while the rest of the chemical 
energy is changed into heat. Under the assumption that the only product is water vapor, 
the low heat value can be used which is used for calculating the heat generated.  
 
4.2.2 RADIATOR 
 The radiator plays an important role in rejecting excessive heat to its ambient 
while it makes the PEM fuel cell operate at its desired constant working temperature for 
any load. The capacity of the radiator should be determined considering the amount of 
heat which is carried by the coolant.  
radiatorcoolant
QQ

=
     
  
  (4- 12) 
Moreover, the required temperature difference, coolant flow rate, and heat transfer 
coefficient of the radiator should be utilized to determine the frontal area of the heat 
exchanger [16]. 
()
, 2
,
[]
coolant p coolant coolant
radiator
coolant stack ambient radiator
mC T
A m
TTh
???
=
??

  
  
  (4- 13) 
The heat transfer coefficient is determined by an empirical value proposed by D.G 
Kroger [16]. It is dependent upon the mass flow rate which is a function of the fan speed. 
34?
?
( )
21
2
()()W mK
1.4495 5.9045 0.1157
radiator fan air
air air
hf fm
mm
?
??
==
=? + ?


  
  
  (4- 14) 
where P
r
 is pressure drop (kPa) and h
rad
 is heat transfer coefficient (W m
-2 
K
-1
).  
 The capacity of the heat exchanger is dependent on the amount of heat which 
should be rejected by the coolant. The characteristic of the heat exchanger are determined 
by the heat transfer area, the temperature difference between inlet, or ambient, and outlet 
temperatures, and the heat coefficient which is controlled by the fan speed. Once the 
mass flow rate is determined, the temperature difference of the coolant can be calculated. 
The temperature difference is a typical design parameter which is usually less than 5?C, 
but not over 10?C. As the coolant mass flow rate is increased, the temperature difference 
decreases. However, this increases the parasitic loss of the system due to the coolant 
pump [17].  
__
__
,
()
rad front st out amb rad
coolant out coolant in
coolant p coolant
ATTh
TT
mC
???
=?
??
??

 
  
  (4- 15) 
where A
rad_front
 denotes the frontal area (m
2
) of the radiator and T
st,out
 is the stack outlet 
coolant temperature which is the same as the radiator inlet coolant temperature.  
 The fan is able to increase the mass flow rate of the air which increases the heat 
transfer coefficient which improves the thermal circuit?s cooling characteristics. However, 
the fan also consumes electric power as it spins and blows air. Moreover, the pressure 
drop across the radiator increases as the mass flow rate increases because the friction 
from the surface of the radiator increases [17].  
35?
?
325.101)396.7512.326( +?=
airr
mP 
  
 
  
  (4- 16) 
where P
r
 is pressure drop (kPa).?
The required power for the fan is calculated as follows using the thermodynamic 
power equation under the assumption that compression is an isotropic process [17]. 
1
1
1
k
k
fan air air amb r
elec fan
PmCpTP
??
???
??
??
??
??=?
??
??

 
 
  
  (4- 17) 
where P
fan
 denotes the electric power (W) of the fan.  
 
4.2.3 WATER RESERVOIR 
 Water used as the coolant needs to be stored in a reservoir for safe circulation and 
minimizing additional heat transfer after cooling down through the heat exchanger. There 
are three major kinds of heat transfer in the reservoir system. One is the heat carried in by 
the coolant, another is transfer between coolant and reservoir surface through conduction 
and the other is transfer between reservoir surface and the ambient air through convection 
[17]. 
( ) ()[]
ambirir
p
rescc
resres
p
resres
TTAhTTCpm
Cpm
t
TT ??+??
?
?
?=
,,

  
  (4- 18) 
where T
res
 is the reservoir outlet coolant temperature at the end of time step (K), T
res
p
 is 
temperature of the reservoir at the beginning of time step (K), ?t is time interval (sec), 
m
res
 is equivalent mass of the coolant in the reservoir (kg), Cp
res
 is specific heat (J kg
-1
 K
-
1
), h is heat transfer coefficient of the reservoir, and T
r,i
 is temperature of coolant into 
reservoir (K). 
36?
?
5 ANALYSIS OF A TWO-CELL STACK 
The characteristic of the stack should be analyzed for better operating condition 
and geometry parameter before the stack is designed and operated. Therefore, the two-
stack is analyzed under the varying operating condition and various geometry parameter.  
?
5.1.1.1 COMPARISON BETWEEN SIMULATIONS AND 
EXPERIMENTS FOR A TWO-CELL STACK 
Using the model, a two-cell stack was simulated to analyze static and dynamic 
behaviors. The I-V characteristics calculated for the stack and the individual cells were 
compared with the experimental data collected [6] at a test station installed at Auburn 
University. Operating temperature was 333.15 K. Stoichiometric ratio was 1.2 and 3 for 
each anode and cathode side. Cathode side is humidified as 100% while anode is supplied 
with dry hydrogen gas. 
Figure 14 shows comparison between the simulated and experimental results [6] 
of the voltage across the stack. The simulated and experimental stack voltages matched 
well, but the deviation of the two cell voltages in the experimental results was larger than 
that of the simulated results. Differences in the two cells could be caused by non-uniform 
characteristics of the individual cells such as reactant distributions and the properties of 
the GDL and the membrane. Details describing the effects of various factors on the cell 
voltages are described in Park, SK and Choe, SY [6][16][27]. 
37?
?
 
0.00 0.05 0.10 0.15 0.20 0.25 0.30 0.35 0.40 0.45 0.50
1.2
1.3
1.4
1.5
1.6
1.7
1.8
 
 
Stack vol
t
age (V)
Current density (A/cm
2
)
 Stack voltage simulated
 Experimental result
 
Figure 14 Comparison of the simulated and experimental two-cell stack and individual 
cell voltages  
 
5.1.1.2 COMPARISON BETWEEN A SINGLE PHASE AND A TWO 
PHASE MODEL FOR A TWO CELL STACK 
The structure of the two-cell stack was the same as that depicted in the Fig. 4 
except that a separator between the two coolant channels was used to minimize the 
coupled cooling effects of adjunct cells on performance. Operating conditions were as 
follows: the stoichiometric ratio for the anode was 1.2 and that for the cathode was 2.0. 
The gas pressure, supply gas temperature, and relative humidity on both sides were 1.0 
bar, 333.15
 
K and 100%, respectively. The average temperature of coolants at the exit of 
cell 1 and 2 for the single-phase model was maintained at 333.15
 
K by controlling the 
38?
?
coolant inlet temperature, and the coolant flow rate was constant at 0.003 (kg s
-1
). The 
two-phase model used the same coolant temperature and flow rate. 
I-V characteristics and temperature distribution of a two-cell stack are illustrated 
in Figure 15. The voltages between cell 1 and cell 2 are larger which can be seen in 
Figure 15(a) when current load is increased. The difference in the cell voltages between a 
single-phase and two-phase is affected by liquid water that is often neglected in previous 
single-phase modeling. When the density of the load current is higher than 0.4 (A cm
-2
), 
the influence of liquid water reduces the performance of the cells and the voltage 
difference between the cells becomes larger. As a result, the voltage drop of the two-
phase model is larger than that of the single-phase model. However, it should be noted 
that the liquid water did not significantly limit the current density in comparison to the 
results of other analyses because the gradient of water concentration, the pressure and 
temperature drop along the flow channel have not been fully considered in this modeling 
[18]. 
Figure 15? (b) shows the temperature distributions in the stack at the given profile 
of current density. When liquid water was considered, the overpotential in the cells was 
higher than when in single phase, and subsequently the associated heat produced was 
more than that observed in a single-phase model. As a result, the working temperature in 
the cells became higher. In addition, the geometrical distribution of the heat sources in 
the stack was asymmetrical. Accordingly, heat transfer rate from the heat sources to the 
coolants was not identical, which caused an asymmetrical temperature profile in that the 
temperature in cell 1 was higher than that in cell 2. The ohmic overpotential was affected 
by the temperature and water content in the membrane. The temperature in the membrane 
39?
?
of cell 2 was higher than that of cell 1. In addition, the cathode gas channel temperature 
of cell 2 was lower than that of cell 1. The saturation vapor pressure was also lower 
which augments the RH in cell 2 and water content. Consequently, the ohmic 
overpotential in cell 2 became lower than that in cell 1 and determined the deviation of 
the amplitude of the two cell voltages.  
0.4 0.6 0.8 1.0 1.2 1.4
0.0
0.2
0.4
0.6
0.8
1.0
1.2
1.4
1.6
1.8
2.0
 
 
Stack vo
ltage
 
(V)
Current density (A/cm
2
)
 One phase model
 Two phase model
?
(a) Calculated polarization curves of single-phase and two-phase models for a 2-cell stack. 
40?
?
0 2 4 6 8 10 12 14 16 18 20
332
333
334
335
336
337
 One phase model
 Two phase model
 
 
T
e
mper
ature (K)
Cell layer
?
(b) Temperature distribution in a two-cell stack 
Figure 15 I-V curve and temperature of a two cell stack under the ramp input of current 
density.  
?
5.1.1.3 PARAMETRIC STUDY FOR OPERATING CONDITION 
?
5.1.1.3.1  TEMPERATURE 
Temperature is a parameter that can be set prior to operation of a PEM fuel cell. 
Generally, higher temperatures can make higher cell voltages. However, each optimized 
temperature could be different and should meet each different PEM fuel cell?s system 
conditions. 
 The inner temperatures of a PEM fuel cell are different across each of the layers 
which constitute the fuel cell stack. Therefore, a temperature measurement is dependent 
41?
?
on the particular location where the measurement is taken. For example, the temperatures 
of the inside and outside are different. The temperatures of the anode catalyst layer and 
cathode catalyst layer are also different. There are some different ways to estimate the 
operating temperature of the fuel cell system. The first method is by measuring the 
temperature of the surface of the PEM fuel cell. The next method is by measuring the 
temperature of the air outlet. The last one is by measuring the temperature of the stack?s 
coolant outlet. All of these methods are an estimation of the fuel cell temperature but are 
not exactly same because the temperature differences means that there is heat transfer 
inward or outward in the PEM fuel cell system. Consequently, determining the 
appropriate spot for measuring temperature should be done so by considering the 
material?s heat conductivity, reactant gas mass flow rate, and coolant mass flow rate. 
 Figure 16 shows the temperature dependent I-V, power and overpotentials? 
characteristic curves which are constituted of cell voltages at different operating 
temperatures of 60, 70 and 80
 o
C. When the air temperature increased, the cell voltages 
tended to increase because of predominantly decreased ohmic overpotential. 
An increase in temperature cause lower open circuit voltage (Eq.3-5). Activation 
overpotential and ohmic overpotential also decrease with increase in temperature 
according to increasing in temperature (Eq. 3-6 and 8). The water vapor pressure in the 
GDL increased at an elevated air temperature, and there is more water vapor because of 
the assumption that there is always saturated vapor. Subsequently, the membrane water 
content increased, and the ohmic overpotential decreased.  
42?
?
However, the rate of the decreasing activation overpotential and ohmic 
overpotential compensate that of the decreasing open circuit voltage. Consequently, the 
cell voltage and power in the simulation results increase. 
0.00.20.40.60.81.01.21.41.6
1.12
1.13
1.14
1.15
1.16
1.17
1.18
1.19
1.20
 
 
Op
en
-c
ir
c
u
it
 v
o
lt
ag
e
,
 V
Current density, A/cm
2
 T = 60
o
C
 T = 70
o
C
 T = 80
o
C
0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6
0.2
0.3
0.4
0.5
0.6
0.7
 
 
A
c
t
i
v
a
t
i
on
 o
v
er
p
o
t
e
nt
ia
l,
 V
Current density, A/cm
2
 T = 60
o
C
 T = 70
o
C
 T = 80
o
C
?
0.00.20.40.60.81.01.21.41.6
0.0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1.0
 
 
C
o
n
c
en
tratio
n ove
r
po
tenti
a
l, V
Current density, A/cm
2
 T = 60
o
C
 T = 70
o
C
 T = 80
o
C
0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6
0.00
0.02
0.04
0.06
0.08
0.10
 
 
O
h
mic o
v
e
r
potential, V
Current density, A/cm
2
 T = 60
o
C
 T = 70
o
C
 T = 80
o
C
?
43?
?
0
20
40
60
80
100
120
140
160
0.00.20.40.60.81.01.21.41.6
0.0
0.2
0.4
0.6
0.8
1.0
 
 
Cel
l
 volt
age,
 V
Current density, A/cm2
 T = 60
o
C
 T = 70
o
C
 T = 80
o
C
 Cel
l
 power,
 W
 
Figure 16 I-V characteristic, power and overpotential for different operating temperature 
 
5.1.1.3.2  STOICHIOMETRIC RATIO 
The reactant gas supplied should be greater than or equal to the amount of 
reactant gas required in order to perform desirably. Usually, the required reactant gas, air, 
and hydrogen can be calculated using Faraday?s Law. 
22
22
2
4
H H
OO
I
mM
F
I
mM
F
=
=


?????????
(5- 1)?
where I is current (A), and F is Faraday constant(C mol
-1
). 
 Figure 17 shows the stoichiometric ratio dependent I-V, power, and overpotentials? 
characteristic curves with different stoichiometric ratio of 2, 3 and 4. When the air 
stoichiometric flow ratio increased, the water was easily removed and the water content 
in the membrane decreased. As a result, the ohmic overpotential increased. At the same 
44?
?
time, the oxygen concentration near the catalyst increased and as a result, the activation 
overpotential decreased. However, the rate of the decreasing activation overpotential and 
concentration overpotential is higher than that of the increasing ohmic overpotential, and 
subsequently, the cell voltages increased. 
0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6
1.15
1.16
1.17
1.18
1.19
 
 
O
p
en
-c
irc
u
it
 v
o
lt
a
g
e
,
 
V
Current density, A/cm2
 Stoich.= 2
 Stoich.= 3
 Stoich.= 4
0.00.20.40.60.81.01.21.41.6
0.2
0.3
0.4
0.5
0.6
0.7
 
 
Act
i
va
ti
on
 ov
er
po
ten
t
i
a
l
,
 
V
Current density, A/cm2
 Stoich.= 2
 Stoich.= 3
 Stoich.= 4
 
0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6
0.0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1.0
 
 
Concentration overpoten
t
ial
,
 V
Current density, A/cm2
 Stoich.= 2
 Stoich.= 3
 Stoich.= 4
0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6
0.00
0.02
0.04
0.06
0.08
0.10
 
 
O
h
m
i
c
 o
v
er
po
te
nt
ia
l,
 
V
Current density, A/cm2
 Stoich.= 2
 Stoich.= 3
 Stoich.= 4
?
45?
?
0
20
40
60
80
100
120
140
160
180
0.00.20.40.60.81.01.21.41.61.82.0
0.0
0.2
0.4
0.6
0.8
1.0
 
 
Cel
l
 volt
age,
 V
Current density, A/cm2
 Stoich. = 2
 Stoich. = 3
 Stoich. = 4
 Cel
l
 power,
 W
?
Figure 17 I-V characteristic, power and overpotential for different stoichiometric ratio 
 
5.1.1.3.3  RELATIVE HUMIDITY 
PEM fuel cells have the good characteristics that allow for operation at lower 
temperatures. Thus, water management is one of the most important issues. Maintaining 
high proton conductivity at the electrolyte membrane through appropriate water 
management plays an important role in producing higher cell voltages and increasing 
performance.  
 Figure 18 shows the relative humidity dependent I-V, power, and overpotentials? 
characteristic curves from the simulation?s results at different relative humidity. 
When the relative humidity increased, the cell voltages tended to decrease due to 
the liquid water saturation. The relative humidity of the air inlet affected the 
concentration of oxygen in the gas flow channel as well as water removal. At a fixed air 
pressure composed of the oxygen, nitrogen and vapor, the increase in relative humidity 
46?
?
led to a high pressure of vapor and a low pressure of oxygen. Therefore, the oxygen 
concentration at the GDL/CL interface decreased, and as a result, the activation 
overpotential increased.  
As the relative humidity at the inlet increased, the water content in the membrane 
tended to increase rapidly, and subsequently, the ohmic overpotential decreased. 
Although the activation overpotential increased, the magnitude of the decrease in ohmic 
overpotential was much higher than that of the increase in the activation overpotential. 
Consequently, the cell voltages decreased. 
 
0.00.20.40.60.81.01.21.41.6
1.156
1.158
1.160
1.162
1.164
1.166
1.168
1.170
1.172
1.174
 
 
O
p
en
-circui
t
 vol
t
a
g
e
,
 V
Current density, A/cm
2
 RH = 0.2
 RH = 0.5
 RH = 1.0
0.00.20.40.60.81.01.21.41.6
0.2
0.3
0.4
0.5
0.6
0.7
 
 
Ac
t
i
vat
i
on overpot
ent
i
a
l
,
 V
Current density, A/cm
2
 RH = 0.2
 RH = 0.5
 RH = 1.0
 
0.00.20.40.60.81.01.21.41.6
0.0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
 
 
Concent
rat
i
on over
po
t
e
nt
i
a
l
,
 V
Current density, A/cm
2
 RH = 0.2
 RH = 0.5
 RH = 1.0
0.00.20.40.60.81.01.21.41.6
0.00
0.02
0.04
0.06
0.08
0.10
0.12
 
 
O
h
m
i
c overpotential
,
 V
Current density, A/cm
2
 RH = 0.2
 RH = 0.5
 RH = 1.0
?
47?
?
0
20
40
60
80
100
120
140
160
0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6
0.0
0.2
0.4
0.6
0.8
1.0
 
 
C
e
ll
 v
o
lt
age
,
 
V
Current density, A/cm
2
 RH = 0.2
 RH = 0.5
 RH = 1.0
 Ce
ll p
o
w
e
r
,
 W
?
Figure 18 I-V characteristic, power and overpotential for different relative humidity 
 
5.1.1.3.4  REACTANT GAS PRESSURE 
PEM fuel cells can be operated at low, atmospheric or high, compressed pressure. 
An increase in pressure is able to improve stack performance, but increases parasitic loss 
as well because power from the stack is used for the compressors. Therefore, the total 
benefit from an increase in pressure should be further examined in terms of the total 
system, B.O.P.  
Figure 19 shows the pressure dependent I-V, power, and overpotentials? 
characteristic curves from the simulation?s results at different inlet pressures. When the 
pressure increases, the cell voltages tended to increase due to the oxygen concentration. 
Larger oxygen concentration cause decrease in activation overpotential and concentration 
overpotential. However, higher mass flow rate from the high pressure cause much easier 
water removal and subsequently, ohmic overpotential increases. However, increase in 
48?
?
ohmic overpotential is less than decrease in activation and concentration overpotential. 
Therefore, the cell voltage increased. 
0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6
1.156
1.158
1.160
1.162
1.164
1.166
1.168
1.170
1.172
 
 
O
p
e
n
-
c
ircuit vol
t
age, V
Current density, A/cm2
 P
ca,in
 = 1bar
 P
ca,in
 = 2bar
 P
ca,in
 = 3bar
0.00.20.40.60.81.01.21.41.6
0.2
0.3
0.4
0.5
0.6
0.7
 
 
A
c
t
i
va
t
i
on over
po
t
ent
i
a
l,
 
V
Current density, A/cm2
 P
ca,in
 = 1bar
 P
ca,in
 = 2bar
 P
ca,in
 = 3bar
?
0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6
0.0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1.0
 
 
Con
c
en
t
r
a
t
i
o
n
 o
v
er
po
t
e
n
t
i
a
l,
 
V
Current density, A/cm2
 P
ca,in
 = 1bar
 P
ca,in
 = 2bar
 P
ca,in
 = 3bar
0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6
0.00
0.02
0.04
0.06
0.08
0.10
 
 
Ohm
i
c o
v
er
po
ten
t
i
a
l
,
 
V
Current density, A/cm2
 P
ca,in
 = 1bar
 P
ca,in
 = 2bar
 P
ca,in
 = 3bar
?
0
20
40
60
80
100
120
140
160
0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6
0.0
0.2
0.4
0.6
0.8
1.0
 
 
Ce
ll v
o
ltag
e, 
V
Current density, A/cm2
 P
ca,in
 = 1 bar
 P
ca,in
 = 2 bar
 P
ca,in
 = 3 bar
 
 Cell
 power,
 
W
?
49?
?
Figure 19 I-V characteristic, power and overpotential for different reactant gas pressure 
 
5.1.1.4  OPTIMIZATION OF OPERATING CONDITION 
There are many ways to optimize the design. As seen below, 9 cases are 
suggested for optimization using the Taguchi method. Figure 20 shows different cell 
voltages depending upon each of the 9 cases. Using the model and simulations, the 9th 
case is selected as the best operating condition. The main factor for improving 
performance is temperature and stoichiometric ratio at the cathode. It causes less 
activation overpotential at the catalyst layer and ohmic overpotential at the membrane.  
 
Case 
Temperature 
(
o
C) 
Stoichiometric
ratio 
Relative 
humidity 
Pressure 
(bar) 
Cell power 
(W) 
1 60 1 20 1 76.67 
2 60 2 50 2 98.33 
3 60 3 100 3 110 
4 70 1 50 3 80 
5 70 2 100 1 98.33 
6 70 3 20 2 110 
7 80 1 100 2 80.83 
8 80 2 20 3 100 
9 80 3 50 1 111.67 
Table 2 Three-level orthogonal arrays for 9 different cell operating conditions 
50?
?
0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4
0.0
0.2
0.4
0.6
0.8
 Case1
 Case2
 Case3
 Case4
 Case5
 Case6
 Case7
 Case8
 Case9  
 
V
o
lt
ag
e (
V
)
Current density (A/cm
2
)
?
Figure 20 I-V characteristic for 9 different cell operating conditions 
 
0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4
0
20
40
60
80
100
120
 Cas1
 Case2
 Case3
 Case4
 Case5
 Case6
 Case7
 Case8
 Case9
 
 
Po
we
r (W)
Current density (A/cm
2
)
?
Figure 21 Cell power for 9 different cell operating conditions 
?
51?
?
5.1.1.5  EFFECTS OF GEOMETRY ON THE PERFORMANCE 
In this section, the effects of the porosity and thickness of the GDL, as well as 
membrane thickness, on the performance were analyzed with the model. Figure 22 shows 
the effect of thickness change of the cathode GDL during a multi-step current test of 0.4, 
0.8 and 1.2 (A cm
-2
). 
Figure 22(a) shows effects of the GDL thickness on the cell voltages. When the 
GDL was thicker, the cell voltage would decrease. As current density increased, the cell 
voltage decreased more rapidly due to the effect of the increasing GDL thickness. 
When the GDL became thicker, the liquid water saturation increased because of 
elevated resistivity for the water transport [19]. As a result, the water content in the 
membrane increased as the GDL became thicker because of the increased liquid water 
saturation, as shown in Figure 22(c). Therefore, the ohmic overpotential decreased. 
Likewise, the increasing liquid water saturation in the GDL makes the oxygen mass 
transport harder than before, and the oxygen concentration at the interface between the 
GDL and the CL decreased. Thus, the activation over-voltages increased as shown in 
Figure 22(b). 
The voltage difference between cell 1 and cell 2 was caused by the decreased 
ohmic overpotential in cell 2 which was affected by the increased membrane water 
content. The membrane?s water content was influenced by cell 2?s operating temperature 
which in turn also affects the voltage difference. 
 
52?
?
 
?
Figure 22 Response of a fuel cell to a multi-step current density at different cathode GDL 
thicknesses. (a) Cell voltage, (b) Activation and ohmic overpotential and (c) Membrane 
water content and oxygen concentration near the cathode catalyst layer.  
?
Figure 23 shows the effects of the porosity of the GDL on the cell voltages, 
membrane water content, and oxygen concentration during a multi-step current test of 0.4, 
0.8 and 1.2 (A cm
-2
). Figure 23 (a) shows the effects on the cell voltages. When the 
porosity became larger, the cell voltages increased because of the smaller pressure drop 
and high rate of oxygen transport which implies an increase of effective diffusion 
60 90 120 150 180 210
0.40
0.45
0.50
0.55
0.60
0.65
0.70
0.4
0.6
0.8
1.0
1.2
1.4
 
Ce
l
l
 vo
ltage 
(V)
Time (sec)
 t
GDL
=0.2mm_Cell 1
 t
GDL
=0.2mm_Cell 2
 t
GDL
=0.4mm_Cell 1
 t
GDL
=0.4mm_Cell 1
 Curre
n
t 
d
ens
it
y (A/cm
2
)
60 90 120 150 180 210
0.00
0.05
0.10
0.45
0.50
0.55
0.60
 
 
O
v
erv
o
lt
ag
e (
V
)
Time (sec)
 t
GDL
=0.2mm_Cell 1
 t
GDL
=0.2mm_Cell 2
 t
GDL
=0.4mm_Cell 1
 t
GDL
=0.4mm_Cell 2Activation overvoltage
Ohmic overvoltage
60 90 120 150 180 210
0
2
4
6
8
10
2
3
4
 Ox
ygen
 conce
ntrat
i
o
n
 (mo
l
/m
3
)
Membrane
 water
 con
tent
Time (sec)
 t
GDL
=0.2mm_Cell 1
 t
GDL
=0.2mm_Cell 2
 t
GDL
=0.4mm_Cell 1
 t
GDL
=0.4mm_Cell 2
 
(b)(a)?
(c)
53?
?
coefficient. In particular, the voltages increased at higher currents as the porosity 
increased because of easier oxygen mass transport when the porosity was larger. As a 
result, the oxygen concentration at the interface between the GDL and CL increased and 
the activation overpotential decreased, as shown in Figure 23 (b) and (c). Thus, the cell 
voltages increased. 
 
?
Figure 23 Dynamic responses of the fuel cells to step changes in current density at 
different cathode GDL porosities. (a) Cell voltage, (b) Activation and ohmic overpotential 
and (c) Membrane water content and oxygen concentration near by the cathode catalyst 
layer.  
60 90 120 150 180 210
0.40
0.45
0.50
0.55
0.60
0.65
0.70
0.4
0.6
0.8
1.0
1.2
1.4
 
Cell v
o
lta
g
e
 
(V)
Time (sec)
 GDL
Porosity
=0.5_Cell 1
 GDL
Porosity
=0.5_Cell 2
 GDL
Porosity
=0.7_Cell 1
 GDL
Porosity
=0.7_Cell 2
 C
u
rren
t
 densit
y (A/
c
m
2
)
60 90 120 150 180 210
0.00
0.05
0.10
0.45
0.50
0.55
0.60
 GDL
Porosity
=0.5_Cell 1
 GDL
Porosity
=0.5_Cell 2
 GDL
Porosity
=0.7_Cell 1
 GDL
Porosity
=0.7_Cell 2
 
 
O
v
erv
o
l
t
ag
e
 (V)
Time (sec)
Activation overvoltage
Ohmic overvoltage
60 90 120 150 180 210
0
2
4
6
8
10
1
2
3
4
5
 
Ox
y
gen concent
rati
o
n
 (mol/
m
3
)
Me
m
b
ran
e
 wa
te
r 
cont
ent
Time (sec)
 GDL
Porosity
=0.5_Cell 1
 GDL
Porosity
=0.5_Cell 2
 GDL
Porosity
=0.7_Cell 1
 GDL
Porosity
=0.7_Cell 2
 
(b)(a)?
(c)
54?
?
?
Figure 24 shows the effects of membrane thickness on the cell voltages, membrane 
water content, and oxygen concentration during a multi-step current of 0.4, 0.8 and 1.2 
(A cm
-2
). Figure 24(a) shows the effects of the membrane thickness on the cell voltages. 
When the membrane thickness increased, cell voltage decreased, and it caused more of an 
effect on cell voltage with a higher current density.  
When the membrane became thicker, there was no change in oxygen 
concentration or activation over-voltage. However, the water content became smaller and 
subsequently the ohmic overpotential increased as shown in Figure 24(b) and (c). The 
rate of the voltage difference between cell 1 and 2 increased as current density increased.  
The analysis showed that the thickness of the membrane plays an important role 
in the improvement of a cell?s performance. As the membrane?s thickness decreases, the 
cell?s performance increases.  
 
60 90 120 150 180 210
0.2
0.3
0.4
0.5
0.6
0.7
0.4
0.6
0.8
1.0
1.2
1.4
 
C
e
ll vo
lt
a
ge (
V
)
Time (sec)
 t
MEM
=0.035mm_Cell 1
 t
MEM
=0.035mm_Cell 2
 t
MEM
=0.055mm_Cell 1
 t
MEM
=0.055mm_Cell 2
 Cu
r
ren
t
 d
ens
it
y (
A
/c
m
2
)
60 90 120 150 180 210
0.0
0.1
0.2
0.45
0.50
0.55
0.60
 t
MEM
=0.035mm_Cell 1
 t
MEM
=0.035mm_Cell 2
 t
MEM
=0.055mm_Cell 1
 t
MEM
=0.055mm_Cell 2
 
 
Over
vo
ltag
e (V)
Time (sec)
Activation overvoltage
Ohmic overvoltage
(b)(a)?
55?
?
?
Figure 24 Dynamic responses of the fuel cells during a step change in current density at 
different membrane thicknesses. (a) Cell voltage, (b) Activation and ohmic overpotential 
and (c) Membrane water content and oxygen concentration near the cathode catalyst layer.  
 
5.1.1.6 OPTIMIZATION OF GEOMETRY PARAMETER 
There are many ways to optimize the design. As seen below, 4 cases are 
suggested for optimization using the Taguchi method. Figure 25 shows different cell 
voltages depending upon each of the 4 cases. Using the model and simulations, the 1st 
case is selected as the best geometry condition. The main factor for improving 
performance is less thickness of GDL which causes decrease in water transportation 
resistance. The easier water transportation to the membrane causes less activation 
overpotential at the catalyst layer and ohmic overpotential at the membrane. 
Case GDL (mm) GDL porosity MEM (mm) Power (W) 
1 0.2 0.5 0.035 106
2 0.2 0.7 0.055 102 
3 0.4 0.5 0.035 99
60 90 120 150 180 210
0
2
4
6
8
10
1
2
3
4
5
 O
x
y
g
en
 c
o
n
c
e
n
t
rat
i
on (
m
ol
/m
3
)
Me
mb
ra
ne
 w
a
t
e
r c
o
nte
n
t
Time (sec)
 t
MEM
=0.035mm_Cell 1
 t
MEM
=0.035mm_Cell 2
 t
MEM
=0.055mm_Cell 1
 t
MEM
=0.055mm_Cell 2
 
(c)
56?
?
4 0.4 0.7 0.055 97 
Table 3 Two-level orthogonal arrays for 4 different geometry parameters 
0.0 0.4 0.8 1.2 1.6
0.0
0.2
0.4
0.6
0.8
 Case1
 Case2
 Case3
 Case4
 
 
V
o
lt
ag
e (
V
)
Current density (A/cm
2
)
?
Figure 25 I-V characteristic for 4 different geometry parameters 
0.0 0.4 0.8 1.2 1.6
0
20
40
60
80
100
 Case1
 Case2
 Case3
 Case4
 
 
Pow
e
r (W
)
Current density (A/cm
2
)
?
Figure 26 Cell power for 4 different geometry parameters 
57?
?
6 STACK DESIGN 
Theoretically, the open circuit voltage of a polymer electrolyte membrane fuel 
cell is 1.23V due to the overpotential. However, real experimental voltage is below 1V. If 
a single cell which produces 1 kW of power is needed, it should generate slightly over 
1000A of current. To do this, the active area should be large enough to generate 1000A. 
Also, the cables between the fuel cell and the load need to be thick enough in order carry 
such a high current and to have a very low ohmic resistance. However, the practical 
solution is to use a multi-cell stack where the cells are connected in series. The active 
area of single cell and cross-sectional area of cables can be reduced by increasing the 
number of cells.  
 
6.1 QUADRUPED ROBOT SPECIFICATION 
The specification for the quadruped robot is based on the walking and 
transporting mode. It has two modes and each mode needs same amount of power, about 
6kW with operation of hydraulic device. The 19Nm maximum motor torque is required at 
maximum speed (1.5m/s) and 28Nm maximum torque is required during the transporting 
mode of the robot. Each rotational speed is about 3000 RPM and 2000 RPM.
58?
?
Condition Robot mass Maximum speed Required power 
Standing up 80 kg 0 m/s 
4.67 kW 
44.6Nm@1000RPM 
Walking mode I 80 kg 1.2 m/s 
8.2 kW 
31.5Nm@2500RPM 
Walking mode II 80 kg 1.5 m/s 
9.5 kW  
33.5Nm@2700RPM 
Transporting mode 130 kg 0.7 m/s 
9.5 kW 
50.2Nm@1800RPM 
Table 4 Operating mode of the quadruped robot 
 
6.2 STACK SIZING 
The first step for the design of the stack is to determine the active area and the 
number of cells in the stack. Furthermore, the desired power output, stack voltage, 
efficiency, and volume and weight limitations are also considered as design requirements. 
Therefore, these variables need to be optimized while the stack constantly satisfies one 
major parameter, such as the desired power output. Cell stack performance is affected by 
the unit cell performance described well by the I-V curve. This is the most important part 
for multi-cell stack sizing, and it can be varied and optimized depending on the operating 
conditions, i.e. pressure, temperature, and humidity. 
 Power output of the stack is calculated by the product of voltage and current.  
stack
P VI=?          (6- 1) 
The total stack voltage is the sum of each unit cell voltage because they are 
connected in series. 
59?
?
1
n
stack i
i
VV
=
=
?
          (6- 2) 
Current is calculated by a product of current density and active area. 
cell
I iA=?           (6- 3) 
There are additional variables, such as stack volume and stack size, which are 
determined by the number of cells in the stack and the active area of each cell. Moreover, 
power output, stack voltage, and stack efficiency can also determine other variables such 
as the number of cells or active area. If there are unknown values among power output, 
stack voltage and stack efficiency, then there are infinite options for the stack design. 
 There are factors that make different combinations of numerous cells and active 
areas difficult to build. For example, connecting a lot of unit cells which have small 
active areas are hard to arrange in a line and connect. Connecting unit cells which have 
larger active areas can also cause problems of uneven distribution of reactant gas 
concentration, temperature, and pressure on the surface of the catalyst layer?s active area. 
One of the most important factors for design of the fuel cell stack is the even distribution 
of reactants, concentration, and pressure. A pressure drop which causes uneven 
distribution should also be fully considered and reduced. Therefore, the final design of 
the fuel cell stack should be checked whether it satisfies all of the desired requirements 
and specifications. If not, the number of cells or active area should be changed 
appropriately to meet the demand. 
 The operating voltage can be any value as long as it meets the specification of the 
system. However, the operating voltage should be determined under the consideration of 
stack size and efficiency. For instance, a fuel cell system operating at a low voltage 
60?
?
makes a higher power output. Thus, as the power density increases, the total stack weight 
and volume can decrease in size for the same power output. However, its efficiency does 
not increase but does decrease which means that the intake components, i.e. fuel tank and 
air blower, should be increased to make up for the poor efficiency in the PEM fuel cell 
system. Inversely, the stack efficiency can be increased with a higher operating cell 
voltage as it allows the PEM fuel cell system to minimize fuel consumptions and the size 
of fuel delivery system. Specifically, in the transportation system, optimizing efficiency 
plays a pivotal role in reducing the weight and size of the system.  
As shown in the Figure 27, the operating cell voltage plays an important role in 
deciding the power density and efficiency. The Figure 28 shows the cell efficiency which 
is based on the relationship between lower heating value of a fuel, hydrogen, and the 
output power. Usually, the power density and the stack efficiency have opposite 
characteristics at the cell voltage larger than 0.5V. The power density increases while the 
stack efficiency decreases, and the power density decreases while the stack efficiency 
increases. The higher operating cell voltage also contributes to make high efficiency and 
reduce fuel consumption even though its size can be increased. Figure 27 shows that the 
stack with operating cell voltage at 0.7V has the 100% more size than the stack which is 
operated at 0.6 V. The stack with operating cell voltage at 0.8V needs also about 4 times 
more active area than the stack with operating cell voltage at 0.6V. Therefore, there are 
both of positive and negative effects about increasing the cell voltage and it should be 
determined based on the priority of the design such as the smallest size or the smallest 
cost for electric generation.  
61?
?
0.0 0.2 0.4 0.6 0.8 1.0
0
2
4
6
8
10
12
14
16
18
20
0.00
0.05
0.10
0.15
0.20
0.25
0.30
 Activareaperwatt
 Power 
densit
y
 (
W
/
c
m
2
)
 
Ac
t
i
v
e
 ar
ea (
c
m
2
/W
)
Cell voltage (V)
 Powerdensity
?
Figure 27 Stack size and operating voltage 
0.0 0.2 0.4 0.6 0.8 1.0
-10
0
10
20
30
40
50
60
70
80
 
 
Cell ef
f
i
c
i
ency
Cell voltage (V)
?
Figure 28 Cell efficiency dependent upon nominal cell potential
 
 
Figure 29 shows the cell power, efficiency and current density of the unit cell 
with the condition that temperature is 353K, RH is 100% and stoichiometric ratio is about 
2. The unit cell with the active area, 250 cm
2
, can make the maximum power, 117 W, at 
62?
?
around 0.5 V with 40% of efficiency. However, the unit cell is able to increase cell 
efficiency up to 55% at nominal cell potential of 0.7V while its maximum power 
decreases about 65%. Therefore, with the condition for the optimization that the cell 
efficiency should be over 45% and the stack should be as small as possible, the operating 
cell voltages can be determined to be 0.6 volt. It generates 98 W with 48% efficiency to 
meet both of performance and efficiency specification. Consequently, 100-cell PEM fuel 
cell system which operates at 0.6 volt to have a 10kW operating power with 50% 
efficiency is needed. 
0.0 0.2 0.4 0.6 0.8
0
20
40
60
80
100
120
0
10
20
30
40
50
60
70
80
90
100
0.0
0.2
0.4
0.6
0.8
1.0
1.2
1.4
1.6
 
 
Cell power (
W
)
Cell voltage (V)
 StackPower
Effi
ci
en
cy
 (
%
)
 Efficiency
C
u
rr
en
t d
e
n
si
ty (
A
/cm
2
)
 CurrentDensity
?
Figure 29 Single cell power, efficiency and current density with selected cell voltage
 
 
6.3 20-CELL STACK ANALYSIS FOR CELL EXTENSION 
?
63?
?
The characteristic of the stack should be analyzed before the single cell is 
expanded to multi cells-stack. Therefore, static behavior and dynamic behavior is 
analyzed to examine a characteristic of temperature distribution within the 20-cell stack. 
 
6.3.1 STACK BEHAVIOR 
First of all, the model was used to analyze the temperature difference under 
constant load current density. Figure 31 shows the simulated results of the voltage 
gradient at constant load current density. It is noted that there is small voltage gradient 
through the stack. It also shows that there is some voltage difference between right end 
cell and left end cell due to the effect of temperature difference. However, in real stack, 
there is much more voltage difference between individual cells [20]. 
0.0 0.2 0.4 0.6 0.8 1.0 1.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
 
 
C
e
ll vo
ltage
 
(V)
Current density (A/cm
2
)
 AverageVoltage
 Cell1
 Cell2
 Cell3
 Cell4
 Cell5
 Cell6
 Cell7
 Cell8
 Cell9
 Cell10
 Cell11
 Cell12
 Cell13
 Cell14
 Cell15
 Cell16
 Cell17
 Cell18
 Cell19
 Cell20
?
Figure 30 I-V polarization curve of the 20 cell stack
 
64?
?
0 5 10 15 20
5.965
5.970
5.975
5.980
5.985
5.990
5.995
6.000
-20 0 20 40 60 80 100 120 140 160
310
312
 Voltage
Cell voltage 
(10
1
V)
Cell number
 Cell Layer
 Te
mp
er
atu
r
e
 (K)
 Temperature
?
Figure 31 Voltage distribution under the constant current density, 0.6 (A cm
-2
) at 2sec
 
 
Figure 32 shows the temperature gradient in the stack at the given current density. 
The peak temperature is observed at the catalyst layer because heat is released at the 
membrane due to the ohmic resistance and at the catalyst layer during chemical reaction 
with the entropy increase. In contrast, the energy loss from the entropy change which 
causes heat generation during the reaction at the anode catalyst layer is negligible. 
Moreover, each catalyst layer exist around membrane at the anodic side and cathode side 
which make temperature difference between right end and left end cell due to 
asymmetrical heat source distribution in the stack [6] [18][25]. 
65?
?
0 1020304050
0.0
0.5
1.0
1.5
 
 
Curr
en
t d
ensi
t
y (A
/
c
m
2
)
Time (sec)
 CurrentDens
-20 0 20 40 60 80 100 120 140 160 180
306
307
308
309
310
311
 
 
 Cell layers
Temp
er
a
t
ur
e
 (
K
)
 2sec
-20 0 20 40 60 80 100 120 140 160 180
315
316
317
318
319
320
321
 
 
Te
mp
er
atu
r
e
Cell Layers
 5sec
-20 0 20 40 60 80 100 120 140 160 180
327.0
327.5
328.0
328.5
329.0
329.5
330.0
330.5
331.0
331.5
 
 
Te
mpe
ra
ture
 (K)
Cell layer
 12sec
?
Figure 32 Temperature distribution under the constant current density 0.6 (A cm
-2
)
 
 
Heat generated from the end cell is transferred to the end plate because of the 
temperature difference between the cell and the endplate which is contact with ambient 
air. Therefore, generated heat from the end cell is not stored in the cell but transferred to 
the end plate assembly. However, the each cathode catalyst layer is located at the right 
side of the membrane assembly in the model. For instance, the cathode catalyst layer of 
the cell 1 is located after 5 layers, coolant channel, gas channel, gas diffusion layer, 
anode catalyst layer and membrane, from the left end plate assembly while there are just 
3 layers between the cathode catalyst layer of the cell 10 and the end plate assembly. It 
causes distribution of heat source to be asymmetrical. Therefore, the amount heat transfer 
66?
?
from the cathode catalyst layer of the cell 10 to the end plate is much larger and it causes 
the right end plate temperature higher than left end plate and temperature gradient in the 
right side relatively smaller than left side. Consequently, overall heat transfer rate at the 
right half side of the stack becomes smaller and it causes to keep the heat and average 
temperature of right half side becomes higher than left side. All in all, heat transfer rate 
from the effect of asymmetrical heat source distribution make asymmetrical temperature 
difference. 
-20 0 20 40 60 80 100 120 140 160 180
305
310
315
320
325
330
335
 
 
Temperature (K)
Cell layer
 2sec
 5sec
 12sec
?
Figure 33 Temperature changes during operation with constant current density
 
 
6.3.2 DYNAMIC BEHAVIOR 
Figure 34
?
shows the current load, the transient behavior of the temperature for the 
each cell which constitutes the stack during the start-up with increasing current load from 
0.1 to 1.5 for 90 seconds. With given operating condition that cell initial temperature is 
298K, relative humidity at the cathode is 70% and anode is 50%, and the stoichiometric 
67?
?
ratio of air is 2 and fuel is 1, it took more than one minute for the cathode catalyst layer 
of the cell 1 (which has lowest temperature) to reach the desired working temperature. 
20 40 60 80 100 120
0.0
0.2
0.4
0.6
0.8
1.0
1.2
1.4
1.6
 
 
C
u
rre
n
t
 de
nsi
t
y (A/
c
m
2
)
Time (sec)
 CurrentDensity
?
0 20 40 60 80 100 120
290
300
310
320
330
340
350
 
 
T
e
m
p
erature (K
)
Time (sec)
 Cell1
 Cell2
 Cell3
 Cell4
 Cell5
 Cell6
 Cell7
 Cell8
 Cell9
 Cell10
 Cell11
 Cell12
 Cell13
 Cell14
 Cell15
 Cell16
 Cell17
 Cell18
 Cell19
 Cell20
?
Figure 34 Temperature changes during the start-up 
68?
?
?
Figure 35 shows the temperature gradient in the stack during the start-up. As 
shown before, the peak temperature is observed at the catalyst layer at which chemical 
reaction occurs and release heat from the entropy increase, while the energy loss at the 
anode catalyst layer with the increase in entropy is negligible. There is also temperature 
difference between right end and left end cell due to asymmetrical heat source 
distribution in the stack which is described in the previous chapter. 
At 60sec, released heat from the cathode catalyst layer increases temperature of 
the cathode catalyst layer and makes peak temperature. The maximum temperature 
difference between the coolant channel and the cathode catalyst area is about 2K at the 
cell 1. It is almost same as the temperature difference between the cell 1 which has lowest 
temperature and the cell 10 which has usually highest temperature. Furthermore, there is 
temperature gradient from the end plate to the 5cells at the both cathode and anode side 
while inner 10 cells showed relatively uniform temperature distribution. 
At 120 sec, there is still peak temperature at the catalyst layer of the unit cell. As 
current density increases, the larger temperature difference between the catalyst layer and 
the coolant channel at the unit cell is observed. The peak temperature has usually more 
than 2K higher temperature than the coolant channel at a current density 1.45 (A cm
-2
). 
Furthermore, temperature difference between cell 1 and cell 20 was about 2K which is 
1K larger than the one at 60. 
69?
?
-20 0 20 40 60 80 100 120 140 160 180
337
338
339
340
341
342
 
 
T
e
m
p
erature (K
)
Cell layer
 60sec
?
-20 0 20 40 60 80 100 120 140 160 180
338
340
342
344
346
348
 
 
T
e
m
per
a
t
ur
e
 (
K
)
Cell layer
 120sec
?
Figure 35 Temperature distribution during the start-up 
 
70?
?
However, the voltage difference from the different temperature is not significant. 
Figure 36 shows the change in cell voltage during the start-up at the current density 0.4, 
0.8 and 1.2 (A cm
-2
). From the simulation, it is shown that excessively high current loads 
cause large temperature gradient from the end plate to the cell 5 at the cathode and anode 
side of the stack and disturb uniform temperature distribution which can improve the 
performance of the stack. Therefore, operating cell voltage should be optimized to realize 
uniform temperature distribution through the stack [6]. Furthermore, asymmetrical 
temperature distribution also should be improved to reduce the difference in heat transfer 
rate at the end cell by using different coolant flow rate or end plate designed for each cell. 
0 5 10 15 20
0.44
0.46
0.48
0.50
0.52
0.54
0.56
0.58
0.60
0.62
0.64
0.66
 
 
Volt
age (V)
Cell Number
 0.4 A/cm
2
 0.8 A/cm
2
 1.2 A/cm
2
?
Figure 36 Voltage changes during the start-up
 
?
Figure 37 shows the average temperature difference percentage of the unit cell in 
the stack during the start-up. From the simulation, there is large temperature gradient at 
both of the end side in percentage terms. There is more than 5% temperature difference 
71?
?
percentage from the 1
st
 cell to the 5
th
 cell and from 17
th
 cell to the 20
th
 cell. Other than the 
cells, the temperature difference between cells can be neglected. 
?
Figure 37 Temperature difference
 
  
6.4 100-CELL STACK 
?
6.4.1 MODEL DESCRIPTION 
Based on those results above, a model for the 100-cell stack was constructed and 
its performance was analyzed. The analysis showed that temperature distribution in the 
stack significantly influenced chemical reactions, phase change between water vapor and 
liquid water, partial pressure and water content in the membranes. However, performance 
of the individual cells was dependent upon the location of cells. The cells located in the 
middle of the stack have not shown any significant gradients of activation as well as 
ohmic over-potentials, partial pressure of the reactants and concentration of water vapor 
0
10
20
30
40
50
60
70
80
90
100
1234567891011121314151617181920
72?
?
and liquid water. Conversely, the outer 5 cells at the cathode and anode side are strongly 
affected by the temperature and showed different performance characteristics from those 
of the middle cells. Maximum voltage difference between the individual cell in the end 
cells and the middle cells was 0.02V at 0.7 (A cm
-2
).  
Cell 
1
Cell 
2
Cell 
3
Cell 
4
Cell 
5
Cell 
96
Cell 
97
Cell 
98
Cell 
99
Cell 
100
?????
?????
Unit 
Cell
?
Figure 38 A schematic diagram of the 100 cell stack 
?
 
Materials 
Thickness
(m) 
Thermal 
conductivity
(W m
-1
 K
-1
)
Density 
(kg m
-3
) 
Specific 
heat 
(J kg
-1
 K
-1
)
End plate G-10 Garolite 0.028575 0.288 1800 1416 
Bus plate 
Stainless steels 
AISI316 
0.001524 13.4 8238 468 
Coolant plate 
AXF-5Q-PYC 
graphite 
0.001524 95 1780 935 
Gas channel 
plate 
TM graphite 0.001524 105 1820 935 
Diffusion media ELAT v3.0 0.0002064 65 558 840 
Catalyst layer Platinum 0.0000035 0.2 387 770 
Membrane 5510 MEA 0.000028 0.21 1967 1100 
Table 5 Geometric parameters and materials for the design 
73?
?
Consequently, a high power stack composed of a high number of cells should be 
constructed considering those end cells. The rest of cells can be assumed to be located in 
the middle of the stack, which can be simply built by a multiplication of the needed 
number of cells with a single cell input and outputs. For building a 100-cell stack, two 
five end cells at each end can be connected to a single cell multiplied with 90, called as a 
stack with a single and multi cells (SSMC). A schematic diagram for the 100-cell stack is 
shown in Figure 38 and Figure 39. Each cell is connected considering temperature 
interaction in series. Reactants flow and chemical reaction has their own characteristic 
and values without interaction with each unit cells, but temperature of the each cell is 
affected by the connection of the cells. 
End plate assembly
Coolant 1
Cell 1
End plate assembly
Coolant 101
Cell 100
Cell 51
Coolant 52
Coolant 51
#
#
Coolant flow
Reactants flow
(Air and hydrogen)
Heat transfer
Coolant out
Coolant in
Reactants in
Reactants out
 
Figure 39 100-cell stack block diagram 
74?
?
?
6.5 DESIGN OF THE AIR PROCESSING SYSTEM 
?
6.5.1 SPECIFICATION OF THE AIR PROCESSING SYSTEM 
There are requirements for each subsystem in order to achieve a desired 
performance. For this fuel cell system, a maximum power of 10.8 kW watts at 233 A and 
46.4V is needed using a 100 fuel cell stack with an active area of 250 cm
2
. Design of the 
blower size and type should be performed carefully because the fuel cell stack receives 
the correct amount of air supply such that the stack can perform at the desired power.  
 First of all, it is needed to estimate how much air should supply the fuel cell stack 
for proper operation.  
13
1
1
233
100 0.0604 ( s )
4 4 96485
0.0604 32 1.93 ( s )
0.0604 22.4 1.35 ( s )
oxygen cell
oxygen
oxygen
I
Nn molm
F
mg
Vl
??
?
?
== =
?
=?=
=?=



    (6- 4) 
Furthermore, the air is composed of 79% oxygen and 21% nitrogen.  
13
1
31 3 3 1
0.79
0.0604 0.0604 0.288 ( s )
0.21
0.288 22.4 6.44 ( s ) 13.7( )
0.00644 ( s ) 1.184 ( ) 7.62 10 ( )
air oxygen nitrogen
air
air
NN N molm
Vlcfm
mmkgm kgs
??
?
????
=+ =+?=
=?= =
=?
 


  (6- 5) 
Pressure drop occurs during operation within the gas channel is the function of the 
friction factor and consumption rate of the reactants. 
75?
?
((), )
ca friction ca consumed
pf mm?= 
     
  (6- 6) 
For the internal flow in the gas channel, pressure drop is defined as follows with 
consideration of friction factor for the developed steady laminar flow in the channel, 
Reynolds number which is relationship between inertial force and viscosity, velocity, 
hydraulic diameter and cross area of gas channel [21]. 
2
56 91
,
2Re
Re
h
h
hD
h
D
L ?v.
?pf wheref
D
VD?
?
==
=
  
   
  (6- 7) 
where ? is the density of the reactants (kg m
-3
) and ? is the viscosity of the reactants (kg 
m
-1
 s
-1
). 
Mass averaged velocity, hydraulic diameter and cross area for the gas channel is 
calculated as shown below [21]. 
4
2( )
c
c
h
c
m
v
A
A
D
P
A ab and P a b
?
=
=
==+

  
 
   
  (6- 8) 
where 
m
is the mass flow of the reactants (kg s
-1
), a and b is each length of the channel 
side (m). 
Once cross area of the gas channel is the square and hydraulic diameter is 
calculated, total pressure drop from the friction can be obtained. 
76?
?
4
28.455
mL
p
a
?
?
??
??
?=
????
??
??

  
 
   
  (6- 9) 
The pressure drop for the reactants which is consumed at the catalyst layer during 
the chemical reaction is as follows 
2
2
2
3
1.93 10
O by reaction
Oconsumed
O
mRT M NRT
p
VV
MN tRT
iA M tRT
VnFV
RT L
Vv
?
? ?
?==
?
?
==
=?

 
 
   
  (6- 10) 
Therefore, total pressure drop is obtained with consideration air and nitrogen 
which constitute air. 
222 2 2
22 2 2
3
4
1.93 10 28.455
air O N O by reaction O along thechannel N
OO N N
ON
pppp p p
mm
RTL L
Vv a
??
??
?
?=?+?=? +? +?
??
????
??
=? + +??????
??
??
????
??

 
  (6- 11) 
 
6.5.2  DESIGN OF THE BLOWER 
A selected blower is required to have a capacity which satisfies both the mass 
flow rate and the pressure drop for the desired maximum operating condition. However, 
different blowers that supply the same mass flow rate provide different pressure rises. 
This means that careful selection or design of a specific blower for this fuel cell system is 
needed.  
There are many kinds of blowers. An axial flow blower and a centrifugal blower 
are the most common. An axial flow blower is proper when much air flow is needed with 
low pressure. However, a high pressure axial flow blower with good efficiency has been 
77?
?
developed recently. According to affinity laws, as impeller diameter increases, air 
volumetric flow rate increases. Also, an increase in fan speed increases static pressure. 
This static pressure increase is relatively greater than the increase of the air volumetric 
flow rate. Therefore, fan engineers usually design a centrifugal fan with a small wheel 
diameter and a high rpm for a pressure blower and an axial flow fan with a large wheel 
diameter and a low rpm for a volume blower. Fan type can be determined considering Ns, 
specific speed. 
 ?
Figure 40 A schematic diagram of blower 
?
6.5.2.1  SPECIFIC SPEED 
The specific speed is constant as long as the fan has geometric similarity and its 
speed remains constant. This is useful for selecting the proper type of pump. For example, 
the specific speed is the speed of any impeller which might be different in size but similar 
Width
Air Inlet
<Side View>
O
u
t
e
r D
i
amet
er
In
n
e
r D
i
amet
er
<Front View>
78?
?
in shape as it delivers 1 cfm and 1 inWC static pressure. It is not related to the scale of 
the blower but by the type and structure of the blower [22]. 
cfm
s
static
m
NN
P
=

      
 
  (6- 12) 
The magnitude of N
s
 varies from 1,000 to 500,000. When it is a small value, it is able to 
boost a large air volume, and when it is a large value, it is able to create a high pressure. 
 
6.5.2.2 TOTAL PRESSURE 
Static pressure is the pressure which affects the boundary parallel with fluid flow. 
It is usually measured at a hole of the surface which is perpendicular to the surface. Total 
pressure is defined as the sum of static pressure and dynamic pressure [22]. 
2
1
2
tsd
d
PPP
Pv?
=+
=
 
      
 
  (6- 13) 
In a fuel cell system, pressure drop occurs from friction within the gas channel, 
thus supplying appropriate air pressure into the gas channel is pivotal for the air supply. 
 
6.5.2.3  AIR BRAKE POWER 
Total air brake power is calculated by considering air flow rate, m , and total 
pressure, P
t
 [22]. 
( ) ( ) 0.056( )
102 60 6120 75 60 4500
tt tt
a
mP mP mP mP
L kW PS PS
?? ??
== == =
 
 
  (6- 14) 
?
79?
?
6.5.2.4  SHAFT HORSEPOWER 
When the fan rotates, the power which is applied to the shaft is shaft horsepower [22], HP. 
( ) 0.09( )
4500
t
s
mP
LPSPS
?
?
==
?

    
 
  (6- 15) 
where ? is efficiency. 
 
6.5.2.5 SHAFT DIAMETER 
The moment applied to the shaft, angular speed and shaft horse power has a 
relationship as below [23]. 
3
2
()
75 75 60 716.2
16
0.73( )
s
TTNTN
L PS
T
dc cm
?
?
??
== =
==
?? ? ?
?
??(6- 16) 
where T is a twisting moment(kg m
-1
), N is a rotation speed (rpm), c is a safety factor and 
? is allowable twisting stress of the shaft. 
The requirement for the air flow rate is 32 cfm of air and 19 in WC at 280 A. 
Therefore, fan horsepower required is 0.096 hp and shaft load will be 0.15 hp, if fan 
efficiency is about 65%. Therefore, the total motor load will be 0.18 hp with a 20% 
mechanical loss. 
 
6.5.2.6  OUTSIDE DIAMETER  
Outside diameter of the fan could be calculated using the equation below with 
specific speed [23]. 
80?
?
2
60
0.15( )
t
u
o
gP
K
Dm
N
?
?
==
?? ? ? ?
?
??(6- 17) 
where P
t
 is total pressure (kg
f
 m
-2
), ? is specific mass of gas (kg
f
 m
-3
), g is gravitational 
acceleration (m s
-2
) and N is speed (rpm). 
 ?
6.5.2.7 INDISE DIAMETER 
From the relationship between outside diameter and inside diameter with given 
specific speed [23],  
0.3
i
o
D
D
=
  
?????
?
??(6- 18) 
Hence, inside diameter is 0.045(m). 
?
6.5.2.8 NUMBER OF BLADES, Z 
The number of blades is calculated with relationship between outside and inside 
diameter as follows [23]. 
2
4sin sin90
8.5 12
1
1.5 1
i
i
o
o
z
D
D
D
D
? ? ?
=?=
??
?
?
??
??
?? ? ?
?
??(6- 19) 
 
6.5.2.9 WIDTH 
2
1
1
4
0.25 0.25 0.045 0.0113( ) 11.3( )
ii
i
DDb
bD m m
?
?=
==?= =
??
?
??(6- 20) 
Generally, 0.3~0.4D
1
 is acceptable, thus b
1
 is14 mm [23]. 
81?
?
 
6.5.3 PARAMETRIC STUDY 
?
6.5.3.1  INFLUENCE OF GEOMETRY PARAMETERS 
Three designs of blower is compared in terms of flow rate and head loss with a 
condition that outside diameter is changed and the shape of blower and inside diameter 
are maintained. The difference in pressure, flow rate and horsepower is recognized in the 
following. 
Specification Design 1 Design 2 Design 3 
N
s
(rpm) 33.34 76.67 153.3 
Efficiency 0.48 0.62 0.65 
Diameter of shaft (m) 1.018 0.742 0.0580 
Wheel diameter (m) 0.2865 0.1450 0.0742 
Inside diameter (m) 0205 0.104 0.0532 
Table 6 Design parameter for the blower 
?
As wheel diameter is increased, overall power is also increased at the same RPM 
while its overall efficiency decreases. Larger wheel diameter is appropriate for large flow 
rate and head loss for air supply system. Therefore, the RPM of the blower can be 
changed while maintaining the same mass flow rate based on the results. Smaller wheel 
82?
?
diameter has less power consumption and slightly better efficiency, but increase in 
efficiency does not increase significantly when the wheel diameter is less than 0.14m. 
??
Figure 41 Brake horse power for a different wheel diameter 
?
Figure 42 Head loss for a different wheel diameter 
?
 
0.5 1 1.5 2
x 10
4
10
-3
10
-2
10
-1
10
0
10
1
RPM
B
r
ak
e 
H
o
r
s
e
 P
o
w
e
r
 (
hp)
Design 1
Design 2
Design 3
10
-2
10
-1
10
0
10
1
0
1000
2000
3000
4000
5000
6000
7000
8000
Flow Rate (CMM)
H
ead Los
s
 (
m
m
A
q)
Design 1
Design 2
Design 3
5000rpm 10000rpm
20000rpm
83?
?
The hub-tip ratio is defined as the inner diameter divided by the outer diameter. 
As inner diameter increases, total pressure can increase much, while flow rate is 
maintained at the same amount. When inner diameter is changed and the outer diameter 
and structure of the blower are fixed, the difference in the head loss is as follows. 
30 40 50 60 70 80 90
6
8
10
12
14
16
18
20
22
He
ad
 Lo
s
s
 (in
 
W
C
)
Hub-tip Ratio (%)
?
Figure 43 Head loss for a same mass flow rate with different hub-tip ratio 
?
6.6 DESIGN OF HEAT EXCHANGER 
Operating conditions 
Max heat generated 22.41 kW 
Air inlet temperature 40
o
C =313.15K null300K 
Water inlet temperature 90
o
C =363.15null365K 
Air inlet speed 10 km h
-1
 = 2.78m s
-1
 
Table 7 Operating condition for the heat exchanger 
84?
?
 
Air inlet temperature is set to be 300K because the quadruped robot system 
traditionally operates under this condition. Furthermore, the PEM fuel cell stack should y 
not be operated over 85
o
C to prevent any degradations. Air flow inlet speed is assumed to 
be 10 (km h
-1
) because the quadruped robot?s speed is 10 (km h
-1
) during normal 
operation. Other properties related to the heat exchanger are as follows [24]: 
Type CF-7.0-5/8 J (External air cooling type) 
Tube outside diameter 16.4 mm =1.64 cm 
Fin pitch 275 / m =2.75 / cm 
Fin thickness 0.0254 cm =0.254mm 
Free flow area / frontal area 0.449 
Heat transfer area / Total volume ?= 269 m
-1 
Fin area / Total area A
f
/A= 0.830 
Air passage hydraulic diameter D
h
=6.68 mm 
Table 8 Parameters for heat exchanger 
?
 
The mass flow rate of water required is calculated using the heat balance 
equation. Therefore, the heat exchanger should be large enough to transfer heat of the 
following amount: 
,
20000
0.475 / sec
4209 10
c
Pc c
Q
mkg
CT
===
?? ?

?
??
?
??(6- 21) 
85?
?
Likewise, the air mass flow rate provided to the heat exchanger and the air 
temperature rise through the heat exchanger can be calculated with the following 
equations: 
1
1.1614 0.3 2.78 0.96( s )
air air fontr air
mAv kg?
?
=? ?= ?? =
?
?
??(6- 22) 
,,
20000
20.7
0.96 1007
air air out air in
air air
Q
TT T K
mCp
?= ? = = =
??


?
?
??(6- 23) 
Thus, the heat exchanger outlet air temperature (K) is equal to 320.7. 
The heat exchanger has 100 cross-tubes, a 0.8 m
2
 frontal area, and a 0.2 m 
thickness. It can be assumed that mass flow is laminar because the Reynolds Number is 
less than 2300. Therefore, the Nusselt number for a tube with constant heat flux is a 
constant, 4.36. 
4.36
cool h
hD
Nu
k
?
==
?? ? ? ?
?
?
??(6- 24) 
where h
cool
 is convective heat coefficient of the coolant (W m
-2
K
-1
), D
h
 is hydraulic 
diameter (m), and k is the conduction heat coefficient (W m
-2
K
-1
). Hence,  
21
4.36 0.677 4.36
179.98( )
0.0164
c
h
k
hWmK
D
? ?
??
== =
?
?
?
??(6- 25) 
Under the assumption that the thickness of the fin is negligible, it can be written 
as: 
()
13.8
1 1 0.83 0.1430
16.4
fin
cool i
air o air
A
AD
AD A
??
??= ?=
??
?? ?
?
?
??(6- 26) 
where 
air
A  is total surface area of the heat exchanger including the fins (m
2
), A
fin 
is total 
surface area of just all of the fins (m
2
), and A
cool 
is total surface area of the inside of the 
tubes (m
2
). Calculating wall conduction resistance by: 
86?
?
ln ln
ln
22
oo
i
iiio
air w
cool
air i
air air
DD
D
DD
AR
Lk A
A D
A A
??
?
?= = = ?
?? ?
?
?
??
??(6- 27) 
62 1
0.0138 16.4
ln 5.025 10 ( )
2 237 13.8
air w
AR mKW
??
?= = ?
?
?
The mass flux, G, is evaluated using the equation below: 
1.1614 2.78
3.2287
1.4544
/
air air front front
air
free free free free front
VA A
m
G
AA A AA
? ? ???
== = = =

?? (6- 28) 
Hence, the Reynolds number is equal to  
7
1.4544 0.00668
Re 525.972
184.6 10
h
GD
?
?
? ?
== =
?
?? ? ?
??(6- 29) 
And from the heat exchanger characteristics [24],  
0.0145
H
j ?
?? ? ?
????
??(6- 30) 
Hence,  
21
2/3 2/3
1.4544 1007
0.0145 0.0145 26.759 ( )
Pr 0.707
air
air
GCp
hWmK
??
? ?
== =
?
??(6- 31) 
To obtain fin efficiency from the efficiency of annular fins of rectangular profile 
[24], 
22
62
2
1
0.000254
0.0285 0.028627
22
0.000254
(0.0285 0.0164) 0.012227
3.105 10
28.5
1.74
16.4
c
c
pc
c
t
rr m
t
L Lm
ALt m
r
r
?
=+= + =
=+= ? + =
== ?
==
?
?
??(6- 32) 
And 
87?
?
1/2
1/ 2
3/2 3/2
6
51.33
0.012227 0.357
0.0263 3.105 10
h
c
p
h
L
kA
?
??
??
==
??
??
??
??
??
?
??(6- 33) 
From the efficiency of annular fins of rectangular profile [24], it follows that ?
f
 ? 
0.89. Hence, overall efficiency becomes: 
() ()91.089.0183.0111
,
=??=??=
f
f
co
A
A
??
??
?
??(6- 34) 
Then, overall heat coefficient is obtatained [18]: 
()
521
1
11 1
11
1.912 10 0.02197 ( )
0.91 51.33
178.98 0.0978
cool W
air O air
cool
cool
air
AR
Uh
A
h
A
mKW
?
? ?
?
=++
??
??
??
=+?+=
?
?
? (6- 35) 
Thus, 
21
45.5( )
c
UWmK
??
=
? ?????
?
? (6- 36) 
To calculate the correction factor F, the following equation is used. 
,,
,,
,,
,,
300 320.7
2.07
355 365
355 365
0.154
300 365
ci co
io
oi ho hi
ho hi
oi
ii ci hi
TT
TT
R
tt T T
TT
tt
P
Tt T T
?
? ?
== = =
?? ?
?
? ?
== = =
?? ?
?? ? ?
?
? (6- 37) 
Hence, from shell and tube correction factor table [24],  
1F?
?? ? ? ?
????
?
? (6- 38) 
Using the Log Mean Temperature Difference (LMTD) methods, 
,, , ,
,, , ,
()()
ln[( ) / ( )]
(365 320.7) (355 300)
49.46( )
ln[(365 320.7) / (355 300)]
hi co ho ci
lm
hi co ho ci
TT TT
T
TT TT
K
???
?=
???
==
?? ? ?
?
? (6- 39) 
88?
?
Therefore, the required length of the heat exchanger is as follows [18]: 
,
20000
1.73( )
45.5 100 0.0164 49.46
lm CF
q
L m
UN DF T??
== =
????
?
?
? (6- 40) 
 
As shown before, the temperature difference is a typical design factor which is a 
proportional inverse to the coolant flow rate required.  
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1.0
1.1
2
4
6
8
10
12
14
16
18
20
0.2
0.4
0.6
0.8
1.0
T
e
m
p
e
r
a
t
u
r
e
 
d
i
ff
e
r
e
n
c
e
, K
R
a
d
i
a
t
o
r
 
f
r
o
n
t
a
l
 
a
r
e
a
,
 
m
2
C
o
o
la
n
t 
fl
o
w
 r
a
te
, 
k
g
/s
e
c
?
Figure 44 The relationship between coolant temperature difference, flow rate and radiator 
frontal area 
 
89?
?
In other words, if a small temperature difference which makes cell temperature 
distribution even is required, the coolant flow rate should be large enough to reject the 
heat generated from the stack as it causes more parasitic loss in the balance of plant. 
However, at the same temperature difference, the coolant flow rate can be minimized by 
increasing the radiator frontal area. The relationship between those parameters is shown 
in the figure below. Moreover, an increase in frontal area with the same coolant flow rate 
allows for a larger temperature difference, i.e. lower coolant outlet temperature, while it 
maintains the same desired operating temperature and decreases parasitic loss. Therefore, 
considering not only parasitic loss, but also the physical size and shape of the heat 
exchanger is pivotal for developing right thermal circuit. 
 
 
?
90?
?
7 100-CELL PEM FUEL CELL SYSTEM 
PEM fuel cell system is made using the 100-cell stack, air supply system and 
thermal circuit system. After 100-cell stack model is developed, air supply system and 
thermal circuit system is connected to estimate power more accurately. 
Fuel?cell
stack
Thermal?
circuit
system
Air?
supply?
system
coolant
m
coolant
m
stack
T
air
m
coolant
T
ca
p
stack
stack
stack
I
V
P
blower
P
pump
P
fan
P
air
T
air
RH
air
SR
?
Figure 45 Schematic diagram of PEM fuel cell system 
?
Auxiliary components consume the power from the stack for operation of a blower, 
coolant pump, radiator fan. The successful operation of B.O.P is important for safe and 
appropriate operation of the stack system.  
91?
?
Static and dynamic performance of the stack with a single and multi cells (SSMC) is 
compared with a model for a stack that is given by simply multiplying a single cell with 
the cell numbers, called as a stack with a single cell (SSC).  
?
7.1.1 STATIC BEHAVIOR 
The model for the stack proposed was used to analyze static and dynamic 
performance. The results of the simulation were compared with the case when the stack 
was simply constructed by a multiplication of the cell numbers with a single cell.  
I-V characteristics of the each cell in the SSMC and SSC are shown in Figure 46. 
There are temperature difference and effect within the stack and they cause voltage 
difference between cells.  
When the current density gets lower, no noticeable differences in the voltages of 
the stack were shown, but the differences in the voltages increased when the current 
density became higher because the more water removal and heat generation rate and 
larger temperature gradient occurs at the high current density. As a result, the maximum 
power of the SSMC became larger than that of the SSC. The voltage difference between 
each stack is about 8V at most. 
92?
?
?
Figure 46 I-V polarization characteristic curve of each cell 
?
0.2 0.4 0.6 0.8 1.0
-10
0
10
20
30
40
50
60
70
80
90
0
2
4
6
8
10
12
 
Stack vol
t
age (V)
Current density (A/cm
2
)
 Stack with middle unit cell
 Stack with single unit cell
Sta
ck po
w
e
r (
k
W)
 
 
Figure 47 Stack power and voltage
 
 
0.0 0.2 0.4 0.6 0.8 1.0
0.2
0.3
0.4
0.5
0.6
0.7
0.8
 
 
 Cell1
 Cell2
 Cell3
 Cell4
 Cell5
 Cell6-95
 Cell96
 Cell97
 Cell98
 Cell99
 Cell100
 Cell1-100
Vol
t
a
ge (
V
)
Current density (A/cm
2
)
93?
?
7.1.2 DYNAMIC BEHAVIOR 
?
7.1.2.1 START-UP MODE 
Temperature behavior is investigated in membranes of individual cells during a 
startup behavior for the 100-cell stack, as shown in Figure 48, where a ramp current 
density was applied. For analysis of dynamic behaviors of the 100-cell stack, the ambient 
temperature is assumed to be 30K and the cooling of the stack is controlled by turning the 
coolant pump on and off, dependent upon a value of the outlet temperature of the coolant 
that was set at 70K. A ramp current with a rate of 0.01 (A cm
-2
 s
-1
) is applied during a 
start-up, while a multi-step current is applied for transient analysis. 
When the current is applied, then the heat inside the stack is generated and 
temperature of the membranes begins to increase, as shown in Figure 48. Because of high 
heat transfer rate of the end cells from the temperature difference with ambient air, 
temperature of the cells located in the middle of the stack was higher than that of the end 
cells. The temperature difference between the SSMC and SSC amounts to around 50K 
during the start-up that took 100 seconds. In addition, the temperature rise of the 
individual cells varies from cell to cell. However, temperature rises in the SSC relatively 
slower than that of the SSMC due to the higher heat transfer rate. The SSC is assumed 
that each cell has a same amount of heat transfer rate and actually single unit cell faces 
the ambient temperature which make large temperature difference which cause higher het 
transfer rate between unit cell and ambient. 
 
94?
?
0 2040608010
300
305
310
315
320
325
330
335
340
345
350
355
360
365
 
 
T
e
m
p
era
t
ure (K)
Time (sec)
 Cell1
 Cell2
 Cell3
 Cell4
 Cell5
 Cell6-95
 Cell96
 Cell97
 Cell98
 Cell99
 Cell100
 Cell1-100
 
Figure 48 Temperature changes during the start-up 
 
During the start-up, temperature gradient increase with the increase in current 
density which causes more heat inside the stack developed by connecting the middle unit 
cell. The interaction of temperature effect between each cell is considered and it makes 
larger temperature gradient at the end cell than the stack made by the single unit cell. The 
SSMC has higher temperature than the SSC because it is connected with cell which is 
hotter than ambient and it causes much less heat transfer and higher temperature at the 
each cell. The temperature difference between the stack with middle unit cell and the 
stack with single cell is about 8K after 40 sec during the start-up. 
95?
?
0 1 2 3 4 5 6 95 96 97 98 99 100 101
303.0
303.5
304.0
304.5
305.0
305.5
306.0
306.5
307.0
 
 
T
e
mperat
ure (K
)
Cell number
 Stack with middle unit cell at 20sec
 Stack with single unit cell at 20sec
 
Figure 49 Comparison of temperature distribution during the start-up 
 
As current density increases, ohmic overpotential has a much greater affect by the 
membrane resistance and heat from the catalyst layer during the higher rate of chemical 
reaction. Moreover, the deviation not only in the each cell between catalyst layer and 
coolant channel, but also in the stack between end cell and middle cells becomes more 
significant due to the larger heat generation in a short time during the start-up. The 
temperature difference between the end cell and middle cell is about 30K after 100sec 
during the start-up. 
 
96?
?
0 1 2 3 4 5 6 95 96 97 98 99 100 101
300
305
310
315
320
325
330
335
340
345
350
355
360
 
 
T
e
m
p
e
r
at
ur
e (
K
)
Cell number
 60sec
 80sec
 100sec
?
Figure 50 Temperature distribution during the start-up 
?
Figure 51 shows a voltage distribution in the SSMC. The middle cells in SSMC 
have a higher voltage than end cells while each cell of the SSC has a same output voltage. 
When the temperature was low, the difference between the SSMC and SSC was smaller 
at the beginning and became larger as temperature rises. The middle cells with higher 
temperature can remove water easily by evaporation and accelerate chemical reaction 
with lower activation energy (Eq.3-6). Therefore, they have less activation over potential, 
ohmic overpotential, concentration overpotential and larger proton conductivity in the 
membrane. In the experiment, each cell has a different cell voltage [20]. In contrast, the 
stack with single cell has a limitation in the reflection of temperature interaction within 
97?
?
multi cell stack and does not explain voltage increases due to the temperature effect in the 
real stack.  
0123456 95 96 97 98 9 10
0.30
0.35
0.40
0.45
0.50
0.55
0.60
0.65
 
 
V
o
lt
age (
V
)
Cell number
Stack with middle unit cell
 20sec
 40sec
 60sec
 80sec
 100sec
Stack with single unit cell
 20sec
 40sec
 60sec
 80sec
 100sec
?
Figure 51 Voltage changes during the start-up 
  
7.1.2.2 TRANSIENT BEHAVIOR 
For transient analysis, responses of oxygen, vapor concentrations, water content 
of membranes in the middle cells and voltages are calculated as a multi-step current 
applied to the stack as shown in the Figure 52.  
98?
?
400 420 440 460 480 500
0.2
0.4
0.6
0.8
1.0
 
 
Cur
r
e
n
t
 de
ns
it
y
 
(
A
/
c
m
2
)
Time(sec)
 Currentdensity
?
Figure 52 Current density applied
 
 
When the current density increases, temperature of the stack is increased as 
shown in Figure 53 due to the higher heat generation at the catalyst layer. The SSMC is 
connected by multi cell facing each other at the both of end side while SSC is composed 
of single unit cell which is in contact with ambient air and has a larger heat transfer rate 
to the endplate. Therefore, heat generated in SSC is transferred easily and the temperature 
effect causes recovery time of the temperature faster than that in the SSMC with a given 
coolant flow rate when the current load suddenly increases. 
99?
?
400 420 440 460 480 500
340
342
344
346
348
350
 
 
T
e
mpe
r
ature (
K
)
Time (sec)
 SSC
 SSMC
?
Figure 53 Transient responses of temperature in the middle cells
 
 
The results below showed different oxygen concentrations caused by the 
difference in the operating temperature, as shown in Figure 54. The oxygen in the multi-
cell tends to be starving by the time an abrupt change of current occurs, which was 
recovered by a controller that helps replenish the oxygen. The stack with middle unit cell 
has a high temperature and it causes less oxygen concentration than the stack with single 
cell in a given pressure and control volume. Furthermore, it causes smaller overshoot of 
oxygen concentration in the SSMC than SSC. 
100?
?
400 420 440 460 480 500
0
2
4
6
0.2
0.4
0.6
0.8
1.0
 
O
x
y
g
e
n
 c
onc
ent
ra
t
i
o
n
 (
m
ol/
m
3
)
Time (sec)
 SSMC
 SSC
 
 Currentdensity
?
Figure 54 Dynamic response of oxygen concentration to the step input 
 
In the stack with middle unit cell, its vapor concentration increases with 
temperature rise because it is assumed that air is humidified with a specific amount of 
ratio and its saturation pressure is increased with increase in the stack temperature. 
Therefore, a rise in saturation pressure makes vapor partial pressure and concentration 
increase.  
101?
?
400 420 440 460 480 500
10.5
11.0
11.5
12.0
12.5
13.0
0.2
0.4
0.6
0.8
1.0
 Cur
r
ent
 dens
it
y
 (
A
/
c
m
2
)
 
Wat
e
r
 v
apor 
conc
ent
r
a
t
i
on (
m
ol/
m
3
)
Time (sec)
 SSMC
 SSC
?
Figure 55 Dynamic response of vapor concentration to the step input 
?
Water contents in the membranes and temperature behaviors are shown in Figure 
56 . When the current density increased, the liquid water saturation is increased due to the 
accelerated chemical reaction. However, there is no significant difference between SSMC 
and SSC. When the current density is increased, the water contents of membrane 
suddenly decrease because of the higher rate of water uptake from the anode to the 
cathode side. The elevated temperature affected the recovery time of the water contents. 
Increasing temperature causes not only higher diffusivity (Eq.3-20) but also higher 
saturation pressure which makes less water activity (Eq.3-23) and water contents (Eq.3-
24). It causes water transport from the membrane to the cathode much harder by electro 
osmotic force (Eq.3-17). Consequently, decrease in water activity makes less water 
transportation out of the membrane and recovery speed of the membrane water contents 
at the membrane in the SSMC slightly faster than that in the SSC. When the current 
102?
?
increases, temperature in the middle cells rises and the maximum difference of the 
temperature was about 2.5K, as shown in the Figure 53. 
400 420 440 460 480 500
0.130
0.135
0.140
0.145
0.150
SSC
 SSMC
 
 
L
i
q
u
i
d
 
w
a
te
r
 sa
tu
r
a
ti
o
n
 (m
o
l
/c
m
3
)
Time (sec)
?
400 420 440 460 480 500
11
12
13
14
15
 
 
Membra
ne water 
co
nt
e
n
t
Time (sec)
 SSC
 SSMC
?
Figure 56 Transient response of liquid water saturation and water content of membranes 
103?
?
in the middle cells 
The ohmic over-potential and the voltage of the stack are shown in Figure 57
?
and 
Figure 58. The activation overpotential and ohmic overpotential of the membranes 
decrease when the temperature rises (Eq. 3-6 and 8) due to the lower activation energy 
and ohmic resistance at the membrane. As a result, the voltage of the SSC that is 
composed of one single cell without temperature effect stack was smaller than that of 
SSMC. 
400 420 440 460 480 500
0.44
0.46
0.48
0.50
0.52
0.54
 
 
Vo
lt
a
ge (
V
)
Time (sec)
 ACT of SSMC
 ACT of SSC
?
( a ) Activation over potential 
104?
?
400 420 440 460 480 500
0.005
0.010
0.015
0.020
0.025
0.030
0.035
0.040
0.045
 
 
Volt
age (V)
Time (sec)
 OHM of SSMC
 OHM of SSC
?
( b ) Ohmic over potential 
Figure 57 Transient responses of activation and ohmic overpotential in the middle cells 
400 420 440 460 480 500
40
45
50
55
60
65
70
75
 
 
S
t
ack vol
t
a
ge (V)
Time (sec)
 SSMC
 SSC
?
Figure 58 Transient responses of PEM fuel cell stack voltage 
105?
?
7.1.3 DYNAMIC BEHAVIOR DURING OPERATION MODE  
After total PEM fuel cell test station is built, the model was simulated under the 
operation mode to verify that the B.O.P system meets the needs for performance. Total 
net power should satisfy specification which needs about 8kW during normal walking 
mode with hydraulic device in the robot. The load profile and corresponding stack 
voltage is as follows. The current load is determined considering operating voltage and 
current, which is around 0.6 V where efficiency is about 50% and maximum current load 
is limited at 1.0A where power output start decreasing and efficiency becomes less than 
40%. As seen in Figure 60, the system generates power output more than 8kW and about 
11kW at most. 
0 40 80 120 160 200 240
50
100
150
200
250
40
50
60
70
 
 
Current (A)
Time (sec)
 Current
 V
o
lt
a
ge (V
)
 Stackvoltage
?
Figure 59 Current load under the operation mode 
106?
?
0 40 80 120 160 200 240
2
4
6
8
10
12
 
P
o
wer (kW)
Time (sec)
 Netpower
 Stackpower
?
Figure 60 Net power and stack power during the operation mode 
 
As seen in Figure 60, there is region where power decreases suddenly because 
increase in current load makes blower motor rotate faster to maintain desired 
stoichiometric air flow ratio at a steady rate. Therefore, blower consumes more power to 
be accelerated with rotational inertia of the blower and motor. Moreover, sudden increase 
in current load causes stack starvation by the time when an abrupt change of current 
occurs as shown in Figure 54 and the shortage causes a sudden decrease in the net power 
as well. Consequently, appropriate margin power and stoichiometric ratio of the oxygen 
should be considered for successful operation. 
107?
?
0 40 80 120 160 200 240
-400
-200
0
200
400
600
800
1000
0
50
100
150
200
 
 
Bl
ower p
o
wer (W
)
Time (sec)
 BLpower
 Coolantpower
 
Co
ol
an
t
 p
o
w
e
r
 (
W
)
?
Figure 61 Parasitic power of air supply system and thermal circuit during operation 
 
The stack temperature is controlled using coolant system including coolant pump 
and radiator fan which start operating when the stack temperature is over the desired 
temperature to minimize parasitic loss. The coolant power suddenly increases at 130 sec 
when the stack temperature is over the desired temperature, 343 K. However, its parasitic 
power gradually decreases as temperature decreases and is close to desired temperature 
because coolant pump is operated by the magnitude of the temperature difference 
between stack and desired temperature. Temperature of each cell in the stack increases 
with current loads until it reaches desired temperature as shown in Figure 62. Each slope 
for increase in temperature is dependent upon the heat generation rate determined by 
magnitude of the current load. 
108?
?
0 40 80 120 160 200 240
290
300
310
320
330
340
350
360
0.0
0.2
0.4
0.6
0.8
1.0
 
 
Temperature (K)
Time (sec)
 Cell1
 Cell2
 Cell3
 Cell4
 Cell5
 Cell6-95
 Cell96
 Cell97
 Cell98
 Cell99
 Cell100
Curre
nt
 density
 (A/cm
2
)
?
Figure 62 Temperature changes during the operation 
?
7.1.4 DESIGN PROGRESS 
Through using different kinds of material and considering the different operating 
condition and geometry, final design of the 100-cell PEM fuel cell system achieved about 
0.44 (W cm
-2
) power density, 32% reduction in weight and 29% reduction in volume for 
use in the quadruped robot as shown in the Figure 63.  
 
 
1
st
 design 2
nd
 design 3
rd
 design 
Materials 
Polymer, copper and 
graphite 
G-10 Garolite, stainless 
steel and AXF-5Q-PYC 
graphite and Tm graphite 
G-10 Garolite, stainless 
steel, AXF-5Q-PYC 
graphite, Tm graphite , 
another graphite and 
different geometry 
Active area (cm
2
) 
250 250 250 
109?
?
Coolant and gas channel 
contact percentage (%) 
Normal Normal Normal 
Stack weight (kg) 28.32 24.01 19.4 
Dimension 
(W x H x L) (cm) 
22 x 25 x 77 22 x 25 x 66 22 x 25 x 53.93 
Table 9 Design progress of the PEM fuel cell system 
?
0.36
0.40
0.44
0.48
0.52
0.56
16
17
18
19
20
21
 
Po
wer Densit
y 
(
W
/
c
m
^
2)
 PowerDensity
 St
ack weight
 
(k
g)
 StackWeight
?
Figure 63 Schematic diagram and design progress of the PEM fuel cell system 
110?
?
8 CONCLUSIONS  
The dynamic model for a 100-cell stack is developed considering temperature and 
two-phase effects in the GDL in which performance was compared with experimental 
data. Simulations were performed to analyze dynamic performance of the stack. First, I-V 
characteristics of the two-cell stack, as well as the individual cells, were compared under 
different operating conditions and geometry. The voltage differences between cell 1 and 
cell 2 were caused by temperature gradients in each layer and subsequently the ohmic 
overpotentials. However, the liquid water did not significantly affect the limit of the 
current density because the gradient of water concentration, the pressure and temperature 
drop along the flow channel have not been fully considered in the modeling. 
Second, the effects of two-phase water on the performance were studied with 
different geometry to analyze the effects of the key design parameters and operating 
conditions including the cell voltages, three kinds of overpotential, water content, oxygen 
concentration and temperature, subsequently the stack is optimized based on the results 
for the two-cell stack. Voltage differences between a single-phase and two-phase model 
were predominantly affected by the activation overpotentials caused by changes in 
oxygen concentration. With a more porous GDL and a thinner GDL and membrane, the 
cell-voltages increased. However, the effects were not significant at a low current density. 
With a higher stoichiometric flow ratio, air temperature, and relative humidity, the cell-
voltages also increased. 
111?
?
Furthermore, each B.O.P components is considered in the modeling to estimate 
more realistic net power from the PEM fuel cell system. Finally, the stack model which is 
appropriate for a high number of a cell was proposed and constructed based on the 20-
cell model analysis and used to analyze static and dynamic characteristics at start-up and 
at varying loads. The 100-cell model takes into account three important properties of a 
PEM fuel cell that includes electrical behavior of the terminal, temperature in the cells 
and characteristics of the two-phase flow. Comparisons between a stack with a single and 
multi cells (SSMC) and a stack with single cell (SSC) show no significant terminal 
behavior, while temperature distribution, water concentrations, and oxygen concentration 
through the cells of SSMC, to name a few, are completely different from that of SSC. 
Design of a stack, components of balance-of-plant, and the associated controls should be 
based on the behavior of the complete stack. 
 Future work will include three phase effect with ice under the subfreezing 
temperature. Start-up time of the PEM fuel cell stack in mobile applications is of great 
importance for commercialization. In particular, when operation of a stack stops and the 
ambient temperature drops under sub-freezing temperature, water residing in the various 
layers of the stack forms ice. As a result, the start-up time at the following start is 
increased by latent hat which absorbs the heat produced, blocking the flow of the 
reactants in the pores and dropping proton conductivity in the membranes. Thus, the 
start-up time directly relates the amount of residing ice, the produced heat that causes the 
ice to thaw, flow characteristic in the porous media determined by ratio of the ice, liquid 
water and water vapor, and proton transport affected by water content in the membrane. 
112?
?
These complex phenomena should be explored to ensure commercialization of the fuel 
cell. 
113?
?
REFERENCES 
[1] M. R. von Spakovsky and B. Olsommer, Fuel cell systems and system modeling 
and analysis perspectives for fuel cell development, Energy Conversion and 
Management 43 (2002) 1249-1257. 
[2] J. Larimie and A. Dicks, Fuel Cell Systems Explained, Wiley, 2003. 
[3] J. Ahn, Design and analysis of air and coolant control, Auburn University, 2007 
[4] H. Meng and C. Y. Wang, Model of two-phase flow and flooding dynamics in 
polymer electrolyte fuel cells, J. Electrochem. Soc. 152 (2005) 1733-1741. 
[5] J. C. Amphlett, R. M. Baumert, R. F. Mann, B. A. Peppley and P. R. Roberge, 
Performance modeling of the Ballard Mark ? solid polymer electrolyte fuel cell. 
?. Mechanistic model development, J. Electrochem. Soc. 142 (1995) 1-8. 
[6] Y. Shan and S. Y. Choe, A high dynamic PEM fuel cell model with temperature 
effects, J. Power Sources 145, 1 (2005) 30-39. 
[7] J. C. Amphlett, R. F. Mann, B. A. Peppley, P. R. Roberge and A. Rodrigues, J. 
Power Sources 61 (1996) 183?188. 
[8] T. E. Springer, T. A. Zawodzinski and S. Gottesfeld, Polymer electrolyte fuel cell 
model, J. Electrochem. Soc. 138 (1991) 2334-2341. 
[9] T. V. Nguyen and R. E. White, A water and heat management model for proton-
exchange-membrane fuel cells, J. Electrochem. Soc. 140 (1993) 2178-2186
 
114?
?
[10] J. M. Correa, An electrochemical-based fuel-cell model suitable for electrical 
engineering automation approach, IEEE 515 (2004) 1103-1113. 
[11] Y. Shan and S.Y. Choe, A segmented model for studying water transport in a 
PEMFC, J. Power Sources 145 (2005) 30?39. 
[12] M. Khandelwal, S.H. Lee and M.M. Mench, One-dimensional thermal model of 
cold-start in a polymer electrolyte fuel cell stack, J. Power Sources 172 (2007) 
816?830. 
[13] J. H. Nam and M. Kaviany, Effective diffusivity and water-saturation distribution 
in single- and two- layer PEMFC diffusion medium, International Journal of Heat 
and Mass Transfer, 46 (2003) 4595-4611. 
[14] J. T. Pukrushpan, A. G. Stefanopoulou and H. Peng, Modeling and control for 
PEM fuel cell stack system, Proceedings of the American Control Conference, 
Anchorage, Alaska, May 8-10, 2002. 
[15] Phoenix Design & Technologies, The PDAT Turbomix. 
http://www.padtinc.com/sales/fuelcell/turbomix/, July 2009 
[16] P. Rodtz, G. Paganelli and L. Guzzella, Optimizing air supply control of a PEM 
fuel cell system, Proceeding of the American Control Conference, Denver, CO, 
June, pp. 2043-9, 2003. 
[17] S. Gurski, Cold Start Effects on Performance and Efficiency for Vehicle Fuel Cell 
System, Master of Science Thesis, Virginia Polytechnic Institute and State 
University, 2002. 
[18] A. Rowe and X. Li, Mathematical modeling of proton exchange membrane fuel 
cells. J. Power Sources 102, 1-2 (1993) 82-96. 
115?
?
[19] L. Wang, A. Husar, T. Zhou and H. Liu, A parametric study of PEM fuel cell 
performances, J. Hydrogen Energy 28 (2003) 1263-1272. 
[20] P. Rodatz, F. Buchi, C. Onder and L. Guzzella, Operational aspects of a large 
PEFC stack under practical conditions, J. Power Sources 128 (2004) 208-217. 
[21] F. M. White, Fluid Mechanics, second ed., McGraw Hill, New York, 1986. 
[22] F. P. Bleier, Fan handbook-selection, application and design, McGraw-Hill, 1997 
[23] H. J. Park, Fluid Machinery Engineering, Il-jin, 2003. 
[24] F. P. Incropera and D. P. Dewitt, Heat and mass transfer, Wiley, 2003. 
[25] H. Ju, H. Meng and C. Y. Wang, A single phase non isothermal model for PEM 
fuel cells. International Journal of Heat and Mass Transfer, 48, 7 (2005) 1303-
1315. 
[26] D. G. Kroger, Radiator Characterization and Optimization, SAE paper 840380, 
1984. 
[27] S. K. Park and S. Y. Choe, Dynamic modeling and analysis of 20-cell PEM fuel 
cell stack considering temperature and two-phase effects, J. Power Sources 179 
(2008) 660-672. 
[28] Q. Yan, H. Toghiani and H. Causey, Steady state and dynamic performance of 
proton exchange membrane fuel cell (PEMFCs) under various operating 
conditions and load changes, J. Power Sources 161 (2006) 492-502.