## 2 Methodology  
Fluid flow and mass transfer in the gas-distribution plates were modeled numerically from a macro-scale perspective. Computational fluid dynamics (CFD) software provided by Multiphysics COMSOL(r) v.3.4 was used to simulate the resulting electrical current distributions in the five flow field designs. The simulation began when a flow field channel design was drawn using the computer-aided design (CAD) function of this software. The geometric model consisted of three regions: the gas-flow channels, a gas-diffusion layer, and a cathode. In the gas diffusion layer, such as a porous carbon cloth, only diffusional processes occurred, and in the cathode, the electrochemical reaction took place. In the following steps according to Figure 1, the assumptions made in the fuel cell modeling were: (1) a three-dimensional (3D) and a single-phase mode; (2) air and water vapor were considered ideal gases; (3) the flow was laminar, steady, and incompressible; (4) the electrical power behavior of the cell was kinetically controlled by an oxygen reduction reaction (ORR); (5) the cathode and gas-diffusion layer were isotropic and homogeneous; (6) the ohmic potential drop in the current collector and gravity were negligible; (7) the process was isothermal; and (8) the gas phase diffuses in the \(x,y,\) and \(z\) axes. The equations employed for our mathematical simulations of the fuel cell are listed in Table 1, and their scalar equations are shown in Table 2. The scalar equations are listed in Table 3. Table 4 details boundary conditions, and the data came from the literature [11; 12], including the molar diffusion volumes [13], the exchange current density [14], and the gas diffusion layer thickness Equation (2). The Navier-Stokes Equation (1) was solved first in order to determine the gas velocity field in the channels. This equation was subjected to boundary conditions, Equations (5)-(9). Besides, the gas species (i.e., \(N_{2}\), \(O_{2}\), \(H_{2}O\)) were taken into consideration. Therefore, a multi-component diffusion (Stefan-Maxwell) Equation (2) was needed to solve the diffusion layer. The Stefan-Maxwell Equation (2) needs to incorporate scalar relations for fluid properties [15; 16]. Furthermore, Equation (2) was solved considering the convective flux Equation (10). Meanwhile, Equation (10) was solved using the velocity in the pore [17] and mass flux vectors [18; 19], and these latter equations were determined by electrode-reaction kinetics [20]. The CFD software employs a finite-element method for solving the partial differential equations (PDE). This method involves dividing the domain of the PDEs of interest into a finite number of linear elements (i.e., the \"mesh\"). The triangular mesh offers better conformity to curvature than the square elements, which is used mostly in areas with right angles, since it reduces the number of mesh elements required. The program allows refining of the mesh, carefully adjusting the type and size of elements on any curvature of the geometry. A mesh consisting of 18,000 computational elements was applied to the domain (see Figure 1D). The iterative computations were terminated once the value of the residues fell to less than \(1\times 10^{-6}\). \begin{table}
\begin{tabular}{c c c} \hline \hline
**Design No.** & **Geometry Type** & **Characteristics of Each Channel Design of Bipolar Plate Geometry** \\ \hline
1 & Serpentine & Consists of serpentine channels with abrupt changes in flow direction. \\
2 & Serpentine & Consists of serpentine channels; the direction of the flow was changed gradually. \\
3 & Serpentine & Consists of serpentine channels with modified dimensions regarding \\  & & Designs 1 and 2 \\
4 & Crisscross & This design is based on continually re-directing the flow in a crisscross pattern. \\
5 & Parallel & This design has parallel and straight channels, and each channel \\  & & must be fed with air at a constant pressure. \\ \hline \hline \end{tabular}
\end{table} 
Table 1: Characteristics of the flow channel designs in the bipolar plates.