Numerical and Theoretical Study of Performance and Mechanical Behavior of PEM-FC Using Innovative Channel Geometrical Configurations

Abstract

Proton exchange membrane fuel cell (PEM-FC) aggregation pressure causes extensive strains in cell segments. The compression of each segment takes place through the cell modeling method. In addition, a very heterogeneous compressive load is produced because of the recurrent channel rib design of the dipole plates, so that while high strains are provided below the rib, the domain continues in its initial uncompressed case under the ducts approximate to it. This leads to significant spatial variations in thermal and electrical connections and contact resistances (both in rib–GDL and membrane–GDL interfaces). Variations in heat, charge, and mass transfer rates within the GDL can affect the performance of the fuel cell (FC) and its lifetime. In this paper, two scenarios are considered to verify the performance and lifetime of the PEM-FC using different innovative channel geometries. The first scenario is conducted by adopting a constant channel height (H = 1 mm) for all the differently shaped channels studied. In contrast, the second scenario is conducted by taking a constant channel cross-sectional area (A = 1 mm²) for all the studied channels. Therefore, a computational fluid dynamics model (CFD) for a PEM fuel cell is formed through the assembly of FC to simulate the pressure variations inside it. The simulation results showed that a triangular cross-section channel provided the uniformity of the pressure distribution, with lower deformations and lower mechanical stresses. The analysis helped gain insights into the physical mechanisms that lead to the FC’s durability and identify important parameters under different conditions. The model shows that it can assume the intracellular pressure configuration toward durability and appearance containing limited experimental data. The results also proved that the better cell voltage occurs in the case of the rectangular channel cross-section, and therefore, higher power from the FC, although its durability is much lower compared to the durability of the triangular channel. The results also showed that the rectangular channel cross-section gave higher cell voltages, and therefore, higher power (0.63 W) from the fuel cell, although its durability is much lower compared to the durability of the triangular channel. Therefore, the triangular channel gives better performance compared to other innovative channels.

1. Introduction

Fuel cell (FC) production is significantly impacted by cell assembly and stacking design [1]. There is a need to adequately press the seam to inhibit the reactants’ penetration and reduce the friction between the layers to maintain any components from the FC stack. The clamping strength remains equal to the force expected for compressing the FC zones appropriately while allowing flow. The properties of the contact interfaces between the segments are affected by the collection pressure, and there will be stack sealing problems such as fuel leakage, inner combustion and unacceptable friction. If an inappropriate or inconsistent collection pressure is used, the porous structure and the gas diffusion passage may be broken due to excess pressure. This can reduce the flow into GDL or destroy MEA in the previous two cases; the cell’s performance may decrease due to the stabilization pressure. Due to variations within the FC elements and stack configuration, as each stack has its own assembly pressure, the numerous primary purposes within stack design and assembly continue to be achieved for proper uniform pressure recognition for thin dimensions and mechanical intensity concerning the electrolyte layer and the electrodes in relation to the gaskets, dipole plates and endplates. To clamp the stack, there is the most common method of using bolts, where the torque factor sometimes needed for bolts and the portions that add to the excellent torque are considered by looking at the optimum stabilization pressure of the characteristics of the FC stack. The optimum torque occurs not only because of the optimal stabilization pressure at the FC layers, but is also further affected in the form and content of the bolt and nuts, the bolt seating, the strings, the layers of the stack and the thickness and the number of layers of the bear loads of the momentum, by clamping the surfaces together to accumulate these materials. The bolt behaves as a strength for resisting the momentum generated through force or external force; the part’s edge performs as a fulcrum [2]. Contact resistance and GDL permeability are regulated through GDL contact and bipolar plate layers’ substantial properties. Since the catalyst and membrane layers are fused, their contact resistance is insufficient. Since the materials are normally non-porous and have similar material properties (great density, related Poisson ratios and Young’s modulus), contact friction between the bipolar plates and certain other layers remains a problem. There are some features regarding the layers of bipolar plate and GDL that make the resistance and permeability of contact more significant than among other layers:

(1)

There are significant variations between the Young’s modulus as well as Poisson ratios (a rigid with smooth materials).

(2)

The layer of GDL remains porous, and due to the drop in pore volume or porosity, the permeability decreases.

(3)

A section of the layer of GDL prevents the flow channels within the bipolar pad, causing more negligible permeability as the compression increases through the GDL [3,4].

To understand the vibration performance and behavior of components of the PEM-FC group, for instance, the membrane, gas diffusion layers, the catalyst layer, and dipole plates required here are to understand the vibration properties so that deformation formation can begin and accelerate. Due to vibration at resonance frequencies, this leads to defects and failure to operate. In addition, the inevitable vibration may lead to an increase in assembly errors during operation of the FC group, especially in cars and their applications, due to the increased vibration and the increase in holes, cracks, and delamination all because of these vibrations, which can lead to poor performance and a decrease in durability [5]. The outcomes of compression strength toward PEM-FCs have been studied in terms of performance. A scanning electron microscope was used to examine the differences in the microstructure of the GDL surface that has been squeezed with the region of the bipolar layer’s rib. The physical properties regarding the compressed region of the GDL are distinct from the uncompressed region just after the construction of the FC. More cracks in the GDL surface of the region are visible on the rib [6]. To study the microstructure of GDL-based carbon fiber materials, a synchrotron X-ray at varying compression degrees for imagery is utilized. The superior clarity of the beam of the synchrotron X-ray and the large spatial resolution regarding the corresponding imaging tool are useful for producing tomographic data that reproduce, with great precision, even the fiber’s most minor morphological information. These data did apply for estimating some of the microstructure’s effects in both gas and liquid transport. Such fibers are twisted near the sides of the flow field channel’s form and turn up with improved concentration rates.

Fiber ends protrude 295 mm through the pipe’s volume at a compression degree of 29 percent, where they theoretically impair the gas flow and shape barriers toward water droplets that are transferred into the flow range tube [7]. The consequences regarding the stiffness of the anisotropic bending of GDLs on the efficiency of FCs are examined, including metallic bipolar plates for the polymer electrolyte membrane. To compare the degree of interference of GDLs, compression was calculated. Variation in the cathode channel’s depth has been found to influence the observed cell activity rather than the anode channel’s depth [8]. A standard carbon paper or fabric GDL remains smooth and the fluid is thus distorted and is imposed within the channel, as the GDL stays squeezed within the plates. Differences in GDL thickness and porosity are due to compression, which influences local transport appearances, as do electrical conductivity, the permeability of gas, and thermal and electrical contact resistances against interfaces with neighboring segments [9]. The flow fields’ geometric parameters may have a tremendous impact on the PEM FC’s overall activity [10]. CFD was used to analyze FC’s performance with varying channel shapes and the diameter of the rib in the direction of the large current operational densities of the PEM. Simulations from the three various channel cross-sections employing related boundary conditions revealed that rib diameter remains one of the most significant factors affecting cell production among the cell’s geometrical parameters. Cell performance declined significantly with reducing rib diameter, even though greater reactant flow valuations and inlet velocities were applied [11]. Zhang et al. [12] carried out an experimental investigation of metal foam materials with more than 90% porosity as alternative flow fields for PEM-FCs, where the effects of the suggested foam metal flow field on the phenomena of transport coupled with the reactions inside the PEM-FC were studied, employing a (3D) non-isothermal with a multiphase model. In the 3D model, the metal foam flow field’s total morphology is taken into consideration after validation with the experimental results. The results show that the metal foam highly enhances the performance of the FC, with the greatest current density distribution with respect to the conventional flow fields of the channel rib. Additionally, the distribution of current and oxygen density within the PEM-FC in the case of the flow field of metal foam shows higher uniformity than in the case of the conventional flow field.

The primary purpose of the present work is to analyze and examine the mechanical performance of the PEM-FC during the assembly process with various engineering components of the channel, but with the same border conditions and to stimulate the division and distribution of pressure inside the cell. A CFD form or design of the PEM-FC has been developed, which is performed during the assembly of the FC (bolt assembling). The PEM-FC simulation design contains the following segments: two bipolar panels, two GDL, two CLs and a membrane. Such a design is applied to investigate and examine the influence of the geometric channel structure toward propagating mechanical stress and pressure within the cell. Design validation is viewed toward the experimental literature’s findings (simulating the FC’s deformation).

2. Mathematical Model of PEM and Numerical Method

2.1. Computational Domain

Developing a mathematical model for a complete cell requires extensive mathematical sources and long simulation periods so that the calculations are limited to this search for a single direct flow channel among land regions. Figure 1 shows the components of the entire mathematical field of the cathode gas flow channels, the anode and the electrode equipment.

Table 1. Dimensions and characteristics of the FC segments.

Property	GDL	Membrane	Bipolar Plate	CL
Material	Carbon Paper	Nafion ®	Carbon Graphite	CCM
Width (mm)	3	3	Channel = 1 m Two ribs = 1 mm each	3
Specific neat (J ${\mathrm{kg}}^{-1}$ ${\mathrm{k}}^{-1}$)	500	1050	750	1050
Poisson’s Ratio	0.25	0.25	0.25	0.25
Expansion Coefficient (${\mathrm{k}}^{-1}$)	$-0.8\times {10}^{-6}$	$123\times {10}^{-6}$	$5\times {10}^{-6}$	$123\times {10}^{-6}$
Density (kg/m3)	400	2000	1800	2000
Channel thickness (mm) with constant cross-sectional area of channel case	0.26	0.24	Thickness of plate: 3 mm	0.03
Channel thickness (mm) with constant height of channel case	0.26	0.24	Height = 1 mm Thickness of plate = 2 mm	0.03

shows the properties and dimensions of the materials for each component [13,14].

The following is assumed for the current simulation:

• Steady-state operating conditions for PEM-FC;
• Constant membrane ionic conductivity;
• The reactive gases permeable through the membrane are neglected as the membrane is impermeable to them;
• Isotropic and homogeneous gas diffusion layer;
• The water entering the cell at the entrance to the cell is neglected;
• The water formed is considered in the vapor phase.

The flow field, momentum and mass transfer equations for the fuel cell channels are:

$$\nabla ·(\mathsf{\rho } \overrightarrow{\mathrm{V}})=0$$

(1)

$$\nabla ·(\mathsf{\rho } \overrightarrow{\mathrm{V}}\times \overrightarrow{\mathrm{V}}-\mathsf{\mu } \nabla \overrightarrow{\mathrm{V}})=-\nabla \left(\mathrm{P}+\frac{2\mathsf{\mu }\nabla \overrightarrow{\mathrm{V}}}{3}\right)+\nabla ·\left(\mathsf{\mu }{(\nabla \overrightarrow{\mathrm{V}})}^{\mathrm{T}}\right)$$

(2)

$$\nabla ·\left(-{\mathsf{\rho } \mathsf{\phi }}_{\mathrm{i}}{{\displaystyle \sum }}_{\mathrm{j}=1}^{\mathrm{j}=\mathrm{N}}{\mathsf{\sigma }}_{\mathrm{ij}}\frac{\mathrm{M}}{{\mathrm{M}}_{\mathrm{j}}}\left[\nabla {\mathsf{\phi }}_{\mathrm{j}}+{\mathsf{\phi }}_{\mathrm{j}}\frac{\nabla \mathrm{M}}{{\mathrm{M}}_{\mathrm{j}}}\right]+{\mathsf{\rho } \mathsf{\phi }}_{\mathrm{i}}· \overrightarrow{\mathrm{V}}\right)=0$$

(3)

2.2. Solid Mechanics Pattern through Assembly

The reason for the large strains in the components of the cell is the compression pressure of the PEM-FCs, where all the components are compressed during the process of cell assembly. This results in a heterogeneous compression load due to the recurring channel rib pattern of the two plates, where the initial state of the area under the channels remains uncompressed, causing spatial differences in GDL porosity and thickness distribution, and in thermal and electrical connections and friction (in each of the interfaces, rib–GDL and membrane–GDL). Cell performance and lifetime are affected by changes in mass, charge, and heat transfer rates through GDL.

Hooke’s law can be adapted to measure the mechanical strain caused in the segments with the stacking and clamping rule [15]:

$$\pi M=\frac{1+v}{E}{\sigma }_{ij}-\frac{v}{E}{\delta }_{ij}{\sigma }_{kk}$$

(4)

where $v$ is the Poisson ratio, $E$ is Young’s modulus, $\sigma$_ij is the stress tensor, $\delta ij$ is the Kronecker delta, and ${\sigma }_{kk}$ = ${\sigma }_{11}$ + ${\sigma }_{22}$ + ${\sigma }_{33}$.

2.3. Compression of PEM-FC’s GDL

It is believed that only the pore volume compresses throughout GDL compression, while its solid volume continues to be fixed. The GDL’s porosity is then decreased and can be estimated as [3]:

$$\epsilon =\frac{{\epsilon }_0-1+e^{vs}}{e^{vs}}$$

(5)

where $\epsilon$₀ represents the initial porosity and the superscript $vs$ shows volumetric strain at each point.

The friction between the gas diffusion layer and the bipolar plate is the primary factor in power loss in PEM-FCs, especially when titanium, molded graphite, or stainless steel are used as the BPP material. Through material properties, stabilization pressure, surface topology, and operating conditions, friction is determined, and as the contact area between the GDL and the dipole plate increases due to the high stabilization compression, friction is reduced. However, increasing the GDL pressure does not increase the flow resistance due to large compression. Furthermore, significant pressure will distort the MEA, resulting in cell leakage and an internal short. As a result, the interaction behavior of the bipolar plate and GDL must be investigated.

The contact resistance within the GDL and bipolar plate, Rcontact (mΩ cm²), is determined by [3]:

$$R_{contact}=2.2163+\frac{3.5306}{P_{contact}}$$

(6)

where $P_{contact}$ is the contact pressure (MPa).

2.4. Computational Grid

The governing equations were discretized and solved using a finite element scheme and solved using a multi-physics computational fluid dynamics (CFD) code. The solutions were subjected to rigorous numerical analyses to guarantee that they were grid independent. It was assumed that when the proportional error was less than 10⁻⁶ within each domain between two consecutive iterations, the solution was deemed to be convergent, and the related set of equations was determined iteratively.

3. Fuel Cell Performance Calculations

The basic structure of a PEM-FC is two electrodes (anode and cathode) separated by a solid membrane that acts as an electrolyte. The membrane is coated on both sides with a thin catalyst layer (Figure 2).

Hydrogen fuel moves within a network of channels into the anode, where it remains divided into protons (at the anode catalyst layer), which pass through the membrane to the cathode, and electrons are captured as an electrical current through an external circuit connecting the two electrodes. The oxidant (air) travels to the cathode through channels where oxygen reacts with electrons within the external circuit and protons passing into the membrane (at the cathode catalyst layer) to produce water. The following chemical reactions take place at the anode and cathode electrodes of a PEM-FC:

Anode reaction:

$$2H_2\to 4H^{+}+4e^{-}$$

(7)

Cathode reaction:

$$O_2+4H^{+}+4e^{-}\to 2H_{2}O$$

(8)

Overall:

$$2H_2+O_2\to 2H_{2}O$$

(9)

While these reactions can have multiple intermediate levels and some (unwanted) side reactions, they correctly reflect the key rules in an FC for the time being. The ultimate reaction is identical to the hydrogen combustion reaction. Combustion is an exothermic reaction, meaning that energy is emitted during the process. Water, DC electricity, and heat are the end products of this process.

Total cell reaction:

$$2H_2+O_2\to 2H_{2}O+electricity+heat$$

(10)

A PEM-FC’s behavior depends largely on the electrochemical reactions that happen with hydrogen and oxygen. It is considered that a system is in thermodynamic equilibrium when three conditions are simultaneously verified: chemical equilibrium, heat balance and mechanical balance.

To explain this situation, we consider the basic thermodynamic functions of temperature, T, pressure, P, volume, V, internal energy, U, and entropy, S. In addition, we must also consider the thermodynamic state functions. These functions’ main characteristic is that their value consists only of the system’s current thermodynamic state, regardless of how the state was achieved. In this sense, several different state functions can be combined to establish another state function, known as Gibbs free energy, G.

Thermodynamic Characteristics and Governing Equations

The enthalpy, H, Gibbs free energy, G, internal energy, U, and entropy, S, are the basic thermodynamic functions. These are first-order homogeneous functions of the system’s components, and they are comprehensive properties of a thermodynamic system. The system’s intensive properties are pressure and temperature, which are zero-order homogeneous functions of the system’s components.

The FC’s electrochemical effects change free energy differences linked by the chemical reaction toward electrical energy in a straight line. The Gibbs free energy difference ($∆G$) within a chemical reaction contains the highest network achievable chemical reaction [16].

$$\Delta G=\Delta H-T\Delta S$$

(11)

where $∆H$ and $∆S$ show the enthalpy change and the entropy change. The potential free energy does equal to the lower enthalpy difference amount $T∆S$, which describes the unavailable energy rising from the system’s entropy difference. The quantity from the heat that does perform with an FC working reversibly is $T∆S$. These functions can be seen in a graphical representation in Figure 3.

The enthalpy (or heat) regarding the chemical reaction signifies the variation among the heats of the development of outcomes and reactants.

$$\Delta H={\left(h_f\right)}_{H_{2}O}-{\left(h_f\right)}_{H_2}-\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$2$}\right.{\left(h_f\right)}_{O_2}$$

(12)

Likewise, as $∆H$, concerning the reaction, signifies the variation within the heats of the production of outcomes and reactants, $∆S$ equals the contrast among entropies regarding products and reactants:

$$\Delta S={\left(s_f\right)}_{H_{2}O}-{\left(s_f\right)}_{H_2}-\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$2$}\right.{\left(s_f\right)}_{O_2}$$

(13)

Generally, electrical work yields the following results in terms of charge and potential:

$$W_{el}=q·E$$

(14)

where $W_{el}$ is the electrical work (J/mole), $q$ is the charge (Coulombs/mole), and $E$ is the potential (E.M.F. of the basic FC) (Volts).

The entire charge conveyed into an FC reaction per mole of H₂ employed is equal to:

$$q=n_e·N_{avg}·q_{el}$$

(15)

Electrical work occurs accordingly:

$$W_{el}=n_e·F·E$$

(16)

As specified earlier, the maximum measure concerning electrical energy produced within the FC resembles Gibbs free energy, $∆G$:

$$W_{el}=-\Delta G$$

(17)

The FC’s theoretical potential is thus:

$$E=\frac{-\Delta G}{n_e·F}$$

(18)

Thus, $∆H$ and $∆S$ denote functions regarding temperature. The specific enthalpy toward any species is displayed through the reaction, which is provided by:

$$h=h_{298.15}+{\displaystyle \underset{298.15}{\overset{T}{\int }}{C}_{p}dT}$$

(19)

Through constant pressure, the particular entropy concerning temperature $T$ exists as follows:

$$s=s_{298.15}+{\displaystyle \underset{298.15}{\overset{T}{\int }}\frac{1}{T}{C}_{p}dT}$$

(20)

where $C_p$ shows the specific heat (J/(mole∙K)), which is also a function of temperature of any gas. It is possible to adopt the following empirical relationship [17]:

$$C_p=a+bT+c{T}^2$$

(21)

where the empirical coefficients a, b, and c, listed in

Table 2. The $C_p$ ’s temperature dependency coefficients.

	a	b	c
${\mathrm{O}}_2$	25.84512	0.012987	−3.9 × 10−6
${\mathrm{H}}_2$	28.91404	−0.0008400	2.01 × 10−6
${\mathrm{H}}_2{\mathrm{O}}_{(\mathrm{gas})}$	30.62644	0.009621	1.18 × 10−6

, are O₂, H₂, and H₂O [2].

Substituting Equation (21) into Equations (18) and (19), and applying to Equations (12) and (15), integrates yields:

$$\Delta H=\Delta H_{298.15}+\Delta a\left(T-298.15\right)+\Delta b\frac{\left(T^2-{298.15}^2\right)}{2}+\Delta c\frac{\left(T^3-{298.15}^3\right)}{3}$$

(22)

$$\Delta S=\Delta S_{298.15}+\Delta a\ln \left(\frac{T}{298.15}\right)+\Delta b\left(T-298.15\right)+\Delta c\left(\frac{T^2-{298.15}^2}{2}\right)$$

(23)

where the differences between the coefficients a, b, and c are $\Delta a$, $\Delta b$, and $\Delta c$, respectively, for reactants and products. These coefficients can be calculated by:

$$\Delta a=a_{H_{2}O}-a_{H_2}-1/2a_{O_2}$$

(24)

$$\Delta b=b_{H_{2}O}-b_{H_2}-1/2b_{O_2}$$

(25)

$$\Delta c=c_{H_{2}O}-c_{H_2}-1/2c_{O_2}$$

(26)

The difference within Gibbs free energy ($G$) can be seen to be during an isothermal rule [18]:

$$dG=vdP$$

(27)

where P is the pressure (Pa), and v is the molar volume (m³/mole).

For an ideal gas:

$$P·v=R·T$$

(28)

Therefore,

$$dG=R·T·\frac{dP}{P}$$

(29)

After integration:

$$G=G_o+R·T·\ln \left(\frac{P}{P_o}\right)$$

(30)

For any chemical reaction:

$$jA+kB\to mC+rD$$

(31)

The difference within Gibbs free energy exists in the difference among outcomes and reactants:

$$\Delta G=m{G}_c+r{G}_D-j{G}_A-k{G}_B$$

(32)

By substituting into Equation (30) [19],

$$\Delta G=\Delta G_o+R·T·\ln \left[\frac{a_C^m·a_D^r}{a_A^j·a_B^k}\right]$$

(33)

where the a are the thermodynamic motion coefficients toward the reacting varieties.

Substituting Equation (18) into Equation (33) gives the relation (the Nernst equation):

$$E=E_o+\frac{R·T}{n_{e}F}·\ln \left[\frac{a_C^m·a_D^r}{a_A^j·a_B^k}\right]$$

(34)

Additionally, the reaction itself is not restricted through the outcome of water concentration on this surface. Therefore, the Nernst equation becomes:

$$E=E_o+\frac{R·T}{n_{e}F}·\ln \left[\frac{\left(P_{H_2}/P_o\right)·{\left(P_{O_2}/P_o\right)}^{0.5}}{\left(P_{H_{2}O}/P_{sat}\right)}\right]$$

(35)

The molar fraction of the species (x_i) is merely the ratio of the partial pressure of the species of interest (P_i) and the total pressure (P_total):

$$x_i=\frac{P_i}{P_{total}}$$

(36)

Substituting Equation (36) in Equation (35), the Nernst Equation for Hydrogen Oxygen/Air FC becomes:

$$E=E_o+\frac{R·T}{n_{e}F}·\ln \left[\frac{\left(x_{H_2}·P_{anode}/P_o\right)·{\left(x_{O_2}·P_{cathode}/P_o\right)}^{0.5}}{\left(x_{H_{2}O}·P_{cathode}/P_{sat}\right)}\right]$$

(37)

The species chosen in the Nernst equation pressure term are only the species which participate in the overall electrochemical reaction. For example, the nitrogen in the cathode is not represented directly, since it does not participate in the electrochemical reaction. Inert species do have an effect on the Nernst voltage, through reduction of the mole fraction of the active species. Additionally, note that the activity values use the pressure of the electrode where the reacting species is located. That is, even though there is water at the anode, it is inert in terms of the hydrogen oxidation reaction and only a participant in the electrochemical reaction at the cathode [20].

The thermodynamic equilibrium condition regarding the electrochemical mode corresponds to the reversible potential obtained from the Nernst equation. Nonetheless, as the current begins to pass within the cell, the cell potential falls below the reversible potential because of the different forms of overpotential, such as activation, ohmic, and concentration overpotential. Additionally, at open circuits, PEM-FCs endure other losses—for instance, internal flows and fuel crossover and the usual electrochemical overpotentials. A standard polarization trajectory of a PEM-FC is displayed in Figure 4. The standard figure of merit for measuring FC efficiency is the polarization curve, reflecting the cell voltage–current relationship. Voltage versus current density, scaled by geometric electrode area, is widely given so that the results can be increased to various cell sizes. Other losses control the curve, which can be divided into four regions. When the outer circuit in Region 1 is disconnected, fuel crossover loss occurs. At low current densities in Region 2, activation loss dominates. The ohmic loss governs Region 3, and the diffusion overpotential is responsible for Region 4, so the polarization curve falls downwards [21].

Valuable output (electrical energy) remains a concern for an FC simply when a current is drawn, though the exact cell potential ($E_{cell}$) is reduced to its equilibrium thermodynamic potential (E) due to irreversible losses. A deviation of the thermodynamic potential happens to correspond with the cell’s electrical output at a current flow. The variation in equilibrium rate that exists is described as the overpotential and has been assigned the symbol ($\eta$). The overpotentials originate prime of activation overpotential (${\eta }_{act}$), ohmic overpotential (${\eta }_{ohmic}$), diffusion (or concentration) overpotential (${\eta }_{diff}$), and crossover of fuel and oxidizer through the electrolyte $({\eta }_x)$.

The potential of the FC can be determined by deducting the total losses of overpotentials of the thermodynamic potential equilibrium:

$$E_{cell}=E-{\eta }_{act}-{\eta }_{ohmic}-{\eta }_{diff}-{\eta }_x$$

(38)

At an open circuit, the current that passed as a result of the electronic leakage equals the ionic exchange, and the following equation can be developed:

$$E_{OCV}=E\times t_i$$

(39)

where $E_{OCV}$ is the open-circuit voltage (V), E is the equilibrium thermodynamic potential (Nernst equation), and t_i is the transference number, which can be calculated by:

$$t_i=\frac{{\sigma }_i}{{\sigma }_i+{\sigma }_e}$$

(40)

Species crossover of reactants through the electrolyte can be a major issue degrading the performance of a PEFC. The 0.2 V loss typically suffered at the open circuit in a PEFC represents about 20% efficiency loss without a useful electron being drawn.

On a semi-empirical basis, the crossover current density ix can be modeled as a base current density that occurs even under open-circuit conditions. In this approach, the current density at the cathode is modeled as (i + i_x). At open-circuit voltage, the external current (i = 0), but i_x is still active.

$${\eta }_x=\frac{RT}{\alpha F}\ln \left(\frac{i_x}{i_o}\right)$$

(41)

The equations for activation polarization in the cathode and anode are shown by:

$${\eta }_{act,c}=\frac{RT}{{\alpha }_{c}F}\ln \left(\frac{i}{i_{o,c}}\right)$$

(42)

$${\eta }_{act,a}=\frac{RT}{{\alpha }_{a}F}\ln \left(\frac{i}{i_{o,a}}\right)$$

(43)

where R is the universal gas constant (8.314 (J/(mole⋅K)), T is the temperature (K), F is Faraday’s constant (F = 96,485 (Coulombs/electron-mole)), $i_{o,c}$ is the cathode exchange current density, $i_{o,a}$ is the anode exchange current density, α_c is the transfer coefficient’s charge within the cathode, and α_a is the transfer coefficient’s charge within the anode.

The so-called Tafel equation can be employed to demonstrate the activation losses simply as:

$${\eta }_{act}=a+b\log i$$

(44)

where

$$a=-2.3\frac{RT}{\alpha F}\log i_o$$

(45)

And

$$b=2.3\frac{RT}{\alpha F}$$

(46)

A Tafel map of $b$ is a plot of the log of current density versus overpotential for a given reaction, and is employed to determine the exchange current density and charge transfer coefficient experimentally. Tafel plots, such as in Figure 5, implement a visual perception regarding the activation polarization of an FC. They are adapted for measuring the exchange current density, provided with the extrapolated intercept at ${\eta }_{act}=0$, a measure regarding the maximum current that can be obtained by negligible polarization, and the transfer coefficient (from the slope).

Alternate current density i_o implies a significant parameter that has a dominating impact on kinetic losses. The current exchange density remains a highly nonlinear function of the kinetic rate fixed of reaction and the local reactant concentration and can be formed including an Arrhenius pattern as:

$$i_o=i_o^{ref}·\exp \left(\frac{-E_{rae}}{RT}\right){\left(\frac{C_{reac\tan t}}{C_{reac\tan t}^{ref}}\right)}^{\gamma }$$

(47)

where $E_{rae}$ is the reaction activation energy (J/mol).

The Butler–Volmer equation, which combines the reactant concentration with the activation overpotential, is used to measure the local current density:

For the cathode side:

$$i_c=i_{o,c}^{}{\left(\frac{C_{O_2}}{C_{O_2}^{ref}}\right)}^{{\gamma }_{O_2}}\left[\exp \left(\frac{{\alpha }_{a}F}{RT}{\eta }_{act,c}\right)-\exp \left(-\frac{{\alpha }_{c}F}{RT}{\eta }_{act,c}\right)\right]$$

(48)

For the anode side:

$$i_a=i_{o,a}^{}{\left(\frac{C_{H_2}}{C_{H_2}^{ref}}\right)}^{{\gamma }_{H_2}}\left[\exp \left(\frac{{\alpha }_{a}F}{RT}{\eta }_{act,a}\right)-\exp \left(-\frac{{\alpha }_{c}F}{RT}{\eta }_{act,a}\right)\right]$$

(49)

where $C_{H_2}$ is the local hydrogen concentration (mole/m³), $C_{H_2}^{ref}$ is the reference hydrogen concentration (mole/m³), $C_{O_2}$ is the local oxygen concentration (mole/m³), $C_{O_2}^{ref}$ is the reference oxygen concentration (mole/m³), $i_{o,a}$ is the anode exchange current density, $i_{o,c}$ is the cathode exchange current density, R is the universal gas constant (8.314 (J/(mole⋅K))), T is the temperature (K), F is Faraday’s constant, α_c is the transfer coefficient’s charge within the cathode, and α_a is the transfer coefficient’s charge within the anode; ${\gamma }^{H_2}$ and ${\gamma }^{O_2}$ are empirically determined concentration parameters.

The ohmic polarization can be represented as:

$${\eta }_{ohmic}=i{\displaystyle \sum _{j=1}^n{r}_j}$$

(50)

To solve ohmic resistance problems, some basic tools are required. First, Ohm’s law can be written as:

$$V=i·A·R$$

(51)

where R is the resistance, measured in units of ohms ($\mathsf{\Omega }=Js/C^2$), and A is the electrode geometric surface area (m²).

The resistance is a function of the geometry regarding the conducting substance. The resistivity r_res ($\mathsf{\Omega }.m$) is an inherent feature concerning a material linked toward the resistance R ($\mathsf{\Omega }$) through the cross-sectional area of ion travel, A (m²), and the linear path length of ion travel (material thickness), L (m).

$$r_{res}=\frac{R·A}{L}\hspace{1em}[\mathsf{\Omega }.m]$$

(52)

For materials with very low ionic resistivity, the inversely related conductivity $\lambda$ is commonly used. By definition,

$$\lambda =\frac{1}{r_{res}}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\left[\frac{1}{\mathsf{\Omega }.m}\right]=\left[\frac{S}{m}\right]$$

(53)

Additionally, Ohm’s law, Equation (51), can be shown as:

$$V=\frac{i·L}{\lambda }$$

(54)

The contact resistance is a function of the material surface state and roughness and the contact pressure between the materials (values of area’s contact resistance are given in units of $\mathsf{\Omega }.m^2$).

$$r_{contact}=\frac{V_{loss}}{i} [\mathsf{\Omega }.m^2]$$

(55)

Concentration polarization occurs when an electrochemical reaction consumes a reactant quickly at the electrode, resulting in concentration gradients. We previously determined that the Nernst equation describes how the potential changes in the electrochemical reaction change the reactants’ partial pressure and this can be expressed as:

$$\Delta V=\frac{RT}{n_{e}F}\ln \left(\frac{C_B}{C_S}\right)$$

(56)

where C_B is the bulk concentration of the reactant (mole/cm³), and C_S is the concentration of the reactant at the surface of the catalyst (mole/cm³).

According to Fick’s Law, the reactant’s flux proportionally changes with the gradient of concentration:

$$N=\frac{A·D·\left(C_B-C_S\right)}{\delta }$$

(57)

where N is the flux of the reactants (mole/s), D is the diffusion coefficient of the reacting species (cm²/s), A is the electrode’s active area (cm²), and $\delta$ is the diffusion distance (cm).

There is a proportional relationship between the current density and the consumption of the reactant per unit area, according to Faraday’s Law. In a steady-state electrochemical reaction, the rate at which the reactant species is consumed is proportional to the diffusion flux:

$$N=\frac{I}{n_{e}F}$$

(58)

By substituting Equation (58) in (57), the following relationship is obtained:

$$i=\frac{n_e·F·D·\left(C_B-C_S\right)}{\delta }$$

(59)

Current density determines the concentration of reactant on the catalyst surface; the higher the current density, the lower the surface concentration. When the absorption rate exceeds the diffusion rate, and the reactant is absorbed faster than it can reach the surface, the surface concentration reaches zero. The limiting current density is the current density at which this occurs. Since there are no reactants on the catalyst surface, an FC can only generate the limiting current as shown in Figure 6. Hence, for C_s = 0, i = i_L, so the limiting current density can be calculated by:

$$i_L=\frac{n_e·F·D·C_B}{\delta }$$

(60)

By combining Equations (56), (59) and (60), the equation of voltage loss by the polarization of concentration is obtained:

$${\eta }_{diff}=\frac{RT}{n_{e}F}\ln \left(1-\frac{i}{i_L}\right)$$

(61)

Each electrode will have this loss (and a separate limiting current density), although the anode loss is typically negligible for a hydrogen feed due to the high concentrations and mass diffusivity, so that a single expression is often used to represent concentration polarization. If the anode’s contribution to the concentration and activation polarization is ignored, little error would result for air–H₂ systems.

The equations for concentration polarization in the cathode and anode are shown by:

$${\eta }_{diff,c}=\frac{RT}{n_{e,c}F}\ln \left(1-\frac{i}{i_{L,c}}\right)$$

(62)

$${\eta }_{diff,a}=\frac{RT}{n_{e,a}F}\ln \left(1-\frac{i}{i_{L,a}}\right)$$

(63)

The cross-sectional area of a geometrical channel influences reactant flow and significantly impacts FC performance. Furthermore, the inlet velocity at the anode and cathode is affected by these geometrical parameters [22]. The inlet velocities of air and hydrogen fuel into the channels are usually described by the following equations:

$$u_{in,c}={\zeta }_c\frac{I}{4F}{A}_{MEA}\frac{1}{x_{O_2,in}}\frac{R{T}_{in,c}}{P_{c,in}}\frac{1}{A_{ch}}$$

(64)

$$u_{in,a}={\zeta }_a\frac{I}{2F}{A}_{MEA}\frac{1}{x_{H_2,in}}\frac{R{T}_{in,a}}{P_{a,in}}\frac{1}{A_{ch}}$$

(65)

where $\xi$ is the stoichiometric flow ratio, A_MEA is area of the MEA (m²), and A_ch is the cross-sectional area of the flow channel (m²).

The parameters listed in

Table 3. Operating and electrode parameters for the base case.

Property	Value
Fuel	Hydrogen
Oxidant	Air
Temperature, T	353.15 (80 °C)
Pressure, P	3 atm
Cathode charge transfer coefficient, ${\alpha }_c$	1
Anode charge transfer coefficient, ${\alpha }_a$	0.5
Current loss, $i_x$	0.002 A/cm2
Reference exchange current density, io	$3\times {10}^{-6}$ A/cm2
Limiting current density, iL	1.4 A/cm2
Internal resistance, rj	0.15 Ohm.cm2

were selected to result in a realistic PEM-FC operation. The analytical equations listed above were modeled using the MATLAB program.

4. Results and Discussion

4.1. Analysis of PEM-FC through Assembly

The temperature and relative humidity of the PEM-FC assembly conditions are 20 °C and 30%, respectively, where they are set so there is no thermal and swelling strain for all parts of the FC and the membrane. In the assembly process, the clamping forces of the bolt and nut are subjected to a particular region of the end plate. A clamping pressure of 1MPa was used during the base case’s assembly. Table 1 displays the material properties of each component. The findings of the CFD report illustrate the research during the assembly phase of the PEM-FC. In all the figures in this section, the PEM-FC’s deformation is the same as seen in the literature. The research is divided into two different parts. The first is to generate five different channels with different cross-sectional area values by simulating five different cross-sections of channels of the same channel height. The second aspect is to model the identical five-channel configurations with the identical cross-section field, resulting in five various channels with varying channel height. The rib width remains the same for both elements. The bipolar configuration is the most common stack configuration. It is important to install the membrane electrode assemblies and the bipolar plates at specific locations. This simulation is employed to investigate and examine the influence of the geometric channel’s structure on the propagation of mechanical pressure and stress within the FC. Through GDL, mass, charge, and thermal transfer rates, we will present the conclusions of each section in the next two sections.

4.1.1. Channel Cross-Sections by an Equivalent Channel Height

Figure 7 shows, during the assembly process, how displacement is divided in the PEM-FC. In a realistic FC, due to the bipolar plate’s channels, this causes the friction pressure on the GDL to differ from the clamping pressure. In addition, it is difficult to predict the friction behavior in the interface without CFD analysis, because of the bipolar plate’s rounded corners. It was seen that the triangular shape channel gives the minimum total displacement compared to the other channels, where, it was observed that the total displacement was about 2.57, 2.53, 2.51, 2.5, and 2.47 μm for rectangular, circle, polygon, trapezoidal, and triangular shaped channels, respectively. Therefore, the triangular shape gives a better performance of the FC due to the smallest deformation during its operational life.

In Figure 8, the development of the pressure division in the FC is illustrated during the assembly process, and it also shows the figures for the pressure distribution (contour plots) and the cell deformation, as these figures show the influence of the repetitive simulation on the pressure division. In Figure 9, it is shown that during the assembly process, we can see the division of pressures in the FC, which was developed where the numbers show the division of pressure from von Mises (contour plots) and the shape of the cell distortion. The ranges of pressure subjected on each channel were about −1.71–8.24, −1.77–8.33, −1.75–8.38, −1.8–8.4, and −1.85–8.45 MPa for rectangular, circular, polygon, trapezoidal, and triangular shaped channels, respectively. Therefore, it is also evident that the FC with a triangular channel is the most durable and resistant to deformation resulting from the stress and pressure distribution during assembly. Therefore, it is considered the best in performance and lifetime as a result of the durability it possesses compared to other studied channels.

4.1.2. Channel Cross-Sections with the identical Cross-Sectional Area

Due to the major spatial variations in the GDL distribution of porosity and thickness, as well as in the contact resistances and thermal/electrical connections (both at the membrane–GDL and rib–GDL interfaces), it was necessary to further investigate the operational life of FC, which is performed by taking a certain uniform cross-sectional area for all five studied channels. Figure 10 shows during the assembly process for total displacement division in the PEM-FC. In a realistic PEM-FC, due to the channels in the dipole plate, this leads to a difference in friction on the GDL from the clamping pressure and it is difficult to predict the contact behavior in the interface without a CDF analysis. It was observed that the total displacement was about 2.61, 2.59, 2.59, 2.58, and 2.58 μm for rectangular, circle, polygon, trapezoidal, and triangular shaped channels, respectively. Therefore, in this case also, the triangular shape gives a better performance of the FC due to having the smallest deformation during its operational life.

Figure 11 presents the division of pressure in the FC that was advanced through the process of assembly. The shapes show the division of contour plots of pressure and the configuration of the FC’s distortion. The shapes clearly show the influence of the shape of the repetitive channel rib on the pressure distribution. In Figure 12, the modified pressure division in the FC through the process of assembly can be observed, and the numbers show the pressure division from von Mises and the distribution. It was shown that the range of pressure subjected to each channel was about −1.7–7.11, −1.74–7.14, −1.72–7.22, −1.77–7.08, and −1.78–7.15 MPa for rectangular, circular, polygon, trapezoidal, and triangular shaped channels, respectively. Therefore, it is also evident that the FC with the triangular channel is the most durable and resistant to deformation resulting from the stress and pressure distribution during assembly. Therefore, it is considered the best in performance and lifetime as a result of the durability it possesses compared to other studied channels.

During the assembly process, it is possible to observe the distribution and evolution of pressure in the FC, as illustrated in Figure 12. The numbers show the distribution of pressure from von Mises (contour shapes) and illustrate the distortion of the FC.

4.2. The Performance of PEM-FC

The curve’s polarization is the most critical feature of an FC in terms of performance. It can be used for FC sizing, tracking, and diagnostics, among other things [23,24,25]. Figure 13 shows how ohmic losses, the losses of activation polarization, and the losses of concentration polarization are subtracted from the potential equilibrium to generate the cell polarization curve. The activation losses’ cathode and anode are grouped for each other; however, due to the oxygen reduction reaction’s slowness, most of the losses occur on the cathode. Figure 14 shows the power diagram of the cell, which explains that the maximum power is 0.63 watt and occurred at a nominal current density of 1.04 A/cm².

The efficiency of FCs with various geometric shapes and identical operating conditions has been investigated. The cross-sections of the five different channels were simulated, with identical conditions in each case. According to the modeling used for these various geometrical channels, it was shown that the rectangular channel cross-section gave better FC voltages, and higher power (0.63 W) from the FC, despite its lower longevity compared to the triangular channel. Therefore, the triangular channel gives better performance compared to other innovative channels, as shown in Figure 15. The rectangular cross-section configuration channel causes a maximum pressure drop through the channel followed by the channel with the polygon cross-sectional area. These configurations of channels have a higher pressure drop, for the corresponding gas flow rate, than that of other cross-sectional configurations. Therefore, it is obvious that the use of the rectangular channel led to an increase in the consumption of hydrogen at the anode, which, in turn, will increase the power produced from the PEM-FC.

5. Conclusions

The design process of the PEM-FC includes the stabilization of pressure, and the material characteristics of each component have a high impact factor on the performance and quality of the PEM-FC. Anchorage pressure is a milestone test for optimum FC performance and life. Pressure that is too high, too low, or too heterogeneous has a negative effect on the efficiency of the FCs. The bipolar stack hot coin is the most commonly used stack. Dipole plates should be fixed in precise positions and membrane electrode groups, unless the error of the assembly appears to be the ideal alignment of nearby dipole plates and changes the deviation in the position of assembly. The aggregation error greatly affects the maximum stress of von Mises in the PEM-FC and the uniformity of the pressure distribution. In order to improve the performance and strength of the PEM-FC, control and the error of the assembly of the dipole plate must be minimized. The main objective of this research is to test the mechanical behavior of a PEM-FC during the assembly process with different engineering configurations for the channel but the same boundary conditions. A model for simulating the pressure distribution within the CFD cell of a PEM-FC has been modified, and the resulting pressure occurs during FC assembly and latch assembly, as the simulation model of PEM-FC involves the following four parts: two GDLs, two dipole plates a membrane, two CLs. The results showed that the channels’ cross-sections have a substantial influence on the uniform pressure distribution and von Mises maximal stress in the PEM-FC. This model is used to investigate and evaluate the channel’s engineering structure in terms of mechanical pressure and stress distribution inside the cell. The results of the simulation model showed that a triangular channel cross-section gave lower deformations and mechanical stresses.

6. Recommendations

The findings show that the CFD solid mechanics model can provide insight and illumination into many of the physical phenomena that lead to experimentally observed FC results. Numerous extensions and enhancements should be considered in order to develop this model further. The following is a list of recommendations for potential work and improvements:

• The model is extended from steady-state to transient.
• Use a CFD technique to simulate the PEM-FCs in a 3D model.
• Develop this model to contain the operation conditions (thermal stresses) of the PEM-FC.