Enhancing the Specific Power of a PEM Fuel Cell Powered UAV with a Novel Bean-Shaped Flow Field

One of the marketing challenges of unmanned aerial vehicles (UAVs) for various applications is enhancing flight durability. Due to the superior characteristics of proton exchange membrane fuel cells (PEMFCs), they have the potential to reach a longer flight time and higher payload. In this regard, a numerical assessment of a UAV air-cooled PEMFC is carried out using a three-dimensional (3-D), multiphase, and non-isothermal model on three flow fields, i.e., unblocked bean-shaped, blocked bean-shaped, and parallel. Then, the results of single-cell modeling are generalized to the PEMFC stack to provide the power of 2.5 kW for a UAV. The obtained results indicate that the strategy of rising air stoichiometry for cooling performs well in the unblocked bean-shaped design, and the maximum temperature along the channel length reaches 331.5 K at the air stoichiometric of 30. Further, it is found that the best performance of a 2.5 kW PEMFC stack is attained by the bean-shaped design without blockage, of which its volume and mass power density are 1.1 kW L−1 and 0.2 kW kg−1, respectively. It is 9.4% lighter and 6.9% more compact than the parallel flow field. Therefore, the unblocked bean-shaped design can be a good option for aerial applications.


Introduction
Unmanned aerial vehicles (UAVs) have attracted much attention for applications such as mapping at spatial and temporal scales, supporting disaster management, precise agriculture and monitoring forest changes, wildlife observation, wireless coverage, public security, land administration, etc. [1]. UAVs are capable of performing some dangerous or challenging tasks with high mobility, protection, and low cost. Nevertheless, it is essential to select a suitable energy source with an optimized energy management system to operate UAVs efficiently [2].
Utilizing lithium-ion batteries in UAVs can supply the required energy for only a limited flight time due to their low energy density. Lithium-ion batteries' energy density is 260 Wh kg −1 , while the energy density of fuel cells is more than 800 Wh kg −1 [3]. In addition to a comparably smaller size, the lack of a need to recharge and continuous operation with instant refuel capability are other advantages of polymer exchange membrane fuel cells (PEMFCs) [4]. Due to the outstanding features of PEMFCs, they are a more promising solution for increasing the flight time capability of UAVs [5,6]. However, considerable developments in their technology and power density are required to make them ideal for aviation applications [7,8]. Fuel cell heat management is one of the suggested approaches to enhance the efficiency and specific power of fuel cells while ensuring reliable and flexible operation. Existing commercial fuel cells with an output power of more than 100 W require a cooling system [9]. using integrated power sources. Herwerth et al. [22] applied a PEMFC in a small UAV to show the capabilities of the endurance of PEMFC for flight. They employed the aircooled PEMFC for thermal management and applied self-humidification at an operating temperature of 60 • C.
Based on the research background, most studies on aircraft propulsion have been performed as experimental methodologies, and the numerical modeling of PEMFC-driven aircraft has not been investigated so far. The goal of this paper is to increase the power density of air-cooled PEMFC in the propulsion system to reduce its weight and enhance the system efficiency using numerical modeling. Through improving the design and the embedded flow field inside BPs, the weight and volume of the PEMFC stack can be potentially reduced. Although the effect of various flow patterns on water management and PEMFC performance have been examined in the literature, their simultaneous effects on weight and volume reduction have not been studied so far. Due to the importance of PEMFC stack weight and size in aviation applications, we numerically compare three flow fields, i.e., parallel, bean-shaped with blockages, and bean-shaped without obstacles, in terms of power density using a three-dimensional (3-D), multiphase, and non-isothermal model. Several goals have been pursued in this paper, which include a proper cooling of the fuel cell using an air-cooled strategy, uniform distribution of temperature, oxygen, and electrical current density, and reaching the maximum amount of mass and volume power density. Figure 1 demonstrates a schematic view of a PEMFC-driven UAV. Since the air stoichiometry is high, the required air for cooling and cathode side is provided simultaneously. The high-pressure hydrogen is stored in the tank and transfers to the anode side of PEMFC after reducing its pressure. The humidifier is removed at the anode side of the fuel cell, because the considered membrane thickness (Nafion ® 212) is thin, and the impact of the electro-osmotic drag against back diffusion is negligible. Unreacted hydrogen fuel is recirculated by a hydrogen pump to the anode side inlet. Additionally, the purge valve precisely controls the operational performance of PEMFC while the stack voltage drops quicker than the estimated voltage.

Assumption
Computational fluid dynamic (CFD) simulation provides conditions for the accurate study of fuel cell behavior. The problem-solving strategy using CFD simulation includes the following steps: determination of domain structure (physical modeling), discretization

Assumption
Computational fluid dynamic (CFD) simulation provides conditions for the accurate study of fuel cell behavior. The problem-solving strategy using CFD simulation includes the following steps: determination of domain structure (physical modeling), discretization of the computational domain and governing equations, specifying the numerical solution approach, extraction of results, and verification of results with available experimental measurements. In this research, non-isothermal, 3D, and two-phase models are applied to the PEMFC simulation. Given the complexity of fuel cell governing equations, the following simplified assumptions are considered,

Governing Equation
The governing equations for PEMFC modeling include mass, species, momentum, energy, and charge equations associated with electrochemical reactions, which are summarized herein.
∇· ερ where µ, P, and S m are viscosity, pressure, and momentum source term. The momentum source term is zero in the gas channel regions because of a very high diffusion rate within the channel and a high permeability coefficient. This term can be described as below for other porous regions [24].
where K p is the permeability coefficient of porous media. The species conservation equation has the form below: where C k is species concentration and D e f f k is an effective diffusion coefficient, which can be obtained by Equation (5) [25] D e f f k where r s is the pore-blocking saturation exponent, s is water saturation, and D 0 k is the reference diffusion coefficient of kth species. The source term of spices equation in the membrane region has the following form.
where I and F are current and Faraday constant. n d is the coefficient of the electro-osmotic drag. It results in the pulling of water molecules, along with protons, to the cathode side, which can be determined by Equation (7) [26] The following are equations to calculate species source terms in CLs [27].
where j an and j cat are the local current density of anode and cathode sides, which can be represented as follows [28]: where r is pore blockage in the current collector area, R is the universal gas constant, j re f is reference current density, T is temperature, ς is specific active surface, α is charge transfer coefficient, and γ is concentration dependence.
The energy conservation equation within the fuel cell can be stated as [29]: where k e f f is an effective thermal conductivity, which represents effective properties in the porous media, and can be calculated using Equation (14) [26].
S T is the source term of energy. This term is due to the irreversibility of electrochemical reactions, phase-change heat, and ohmic heating, which can be expressed as: where η is over-potential, σ is evaporation coefficient, σ m is membrane conductivity, σ s is electrical conductivity, A FG is the surface area of phase-change per unit volume, X sat is maximum water mass ratio in dry air, X H 2 O(g) is the ratio of water in dry air, and ∆h FG is evaporation heat. Based on charge equations, we can calculate delivered protons through the membrane and transmitted electrons through current collectors, GDLs, and CLs. The following equations represent the charge equation for electrons at the anode and cathode electrodes [30].
where ϕ s and j s are the electrical potential of electrodes and the electrical transfer current, respectively. The charge equation for delivering protons at the membrane/CLs interface can be defined as [31]: where ϕ m is membrane ionic potential and σ m is membrane proton conductivity, which can be written by Equation (22) [16].
where λ is membrane water content and depends on water activity (a), which can be estimated using the empirical relation as [32]: where P H 2 O is water vapor pressure and P H 2 O sat is the pressure of water saturation, which depends on temperature [33].
When the PEMFC works at low temperatures (below 100 • C), a part of the produced water at the cathode side converts to liquid water. It can block the porous regions of GDL and prevent reactants gas transfer to the CL, which can deteriorate PEMFC performance. Hence, the investigation of a two-phase model within the fuel cell is essential.
The model of transport and formation of liquid water in the GDLs and CLs regions has the form below [31].
where subscript l indicates liquid water, and r w is condensation rate in the GDLs and CLs, which can be represented as below [24]: where c r is the constant condensation coefficient of water. If the velocities of gas and liquid inside gas channels are supposed as equal, the capillary diffusion within porous regions can be expressed as [31]: where p c is capillary pressure, which can be obtained based on Leveret's function as [34]: where θ c is the contact angle, and σ sur f is the surface tension. The contact resistance between GDLs and BPs is the resistance against the current transfer from BPs to GDLs. The amount of contact resistance between the BP/GDL interface is calculated based on experimental Equation (30) [35].
where P Assembly is assembly pressure, A Contact is the contact surface between the GDL and BP, and A Assembly is the assembly surface. The maximum fuel cell voltage is obtained when the fuel cell operates under reversible thermodynamic conditions. The actual voltage of a PEMFC is less than the theoretical voltage due to reaction drops. Therefore, the cell's operating voltage can be expressed as a deviation from the ideal voltage, where the irreversible voltage subtracts from the reversible voltage at a specific current density. Reversible voltage can be described using Nernst potential as [29]: where ϕ oc is the open-circuit voltage and H is enthalpy. The voltage at the anode side and cathode side can be attained as below [26]:

Boundary Condition
To find unknown parameters of velocity, pressure, the potential of electrodes, temperature, and species concentration, a set of equations, i.e., mass, species, momentum, and energy, along with potential equations, should be solved. For this purpose, we need to specify appropriate boundary conditions to simulate the computational domain. Since the modeling is conducted based on a single-domain method, boundary conditions should only be identified on external surfaces of the computational domain, and there is no need to adjust them in the interface of various zones of the PEMFC. Elimination of internal boundary conditions causes an increase in the model's precision, reducing computational error, and facilitating the computation method. However, utilizing the single-domain method prolongs the computing time. The performed boundary conditions can divide into four categories; inlet, outlet, symmetry, and wall boundary conditions, which are explained herein. At the entrance of both cathode and anode channels, the mass flow inlet boundary condition is employed, and these inlet mass flow rates have the following forms: .
where ζ is the stoichiometric ratio, M is molecular weight, and x is the molar fraction. At the inlet of the channels, the partial water vapor pressure is equal to the saturation pressure at the inlet reactants' humidification temperature. Therefore, the water molar fraction at the inlet of channels is given by Equation (37).
where RH is relative humidity. The molar fraction of the inlet hydrogen and oxygen can be described by Equations (38) and (39), respectively.
The outlet pressure boundary condition is given in the channels' outlet. On the walls, for velocity, the non-slip boundary condition and zero flux conditions for other variables are considered except for the outlet and inlet of the anode and cathode channels. The inlets of the anode and the cathode sides are adjusted counter-current because of thin membrane thickness and removing humidification at the anode side. When air flows along the cathode channel, water concentration at the cathode side is more than the anode side due to water production. On the other hand, the electro-osmotic drag is more at the entrance of the anode side than the channel outlet due to the high availability of hydrogen and further proton transfer to the cathode side. For this reason, the inlet and outlet sides of the anode and the cathode are adjusted in a counter-current configuration to compensate membrane dryness at the entrance of the anode side by water back diffusion from the cathode side. The high back diffusion mechanism is essential to keep high proton conductivity in the membrane. The symmetric boundary condition is utilized in bean-shaped flow field along the flow channel for one-half of the PEMFC. Finally, a boundary condition for the current collector is introduced, which determines the electrical potential equation. The current collector contains the anode and cathode terminals, the electrical potential at the anode terminal is adjusted to zero, and the electrical voltage at the cathode terminal is adjusted to the operating voltage. The considered boundary conditions in both cathode and anode sides of PEMFC are listed in Table A1.

Solution Procedure
The finite volume approach is utilized to solve PEMFC's governing equations with defined boundary conditions using the ANSYS FLUENT fuel cell module. In this module, the equations of potential, electrochemical, and membrane water transfer are added to the Fluent software using User Defined Functions (UDF). It facilitates the prediction of source terms and electrochemical reactions within the PEMFC. The module also incorporates all fuel cell components, including BPs, GDLs, CLs, flow channels, and the membrane, which can be introduced into the software to solve equations. In this program, pressure and velocity fields attain using the SIMPLE algorithm, and a repetition process employs to solve the equations set. This procedure continues to converge the results with an accuracy of 10 −6 . Their conjunction is the main reason to choose the iterative method to solve equations. For two reasons, the conservation equations of PEMFC are strongly linked to each other. The first reason is that the chemical and physical features of the materials vary and depend on the transfer parameters. The second reason is that conservation equations are coupled to each other using the Butler-Volmer equation. For instance, the Butler-Volmer equation is a function of hydrogen and oxygen concentrations; on the other hand, the species concentrations can be identified using the species equation. The Butler-Volmer equation enters as a source term into the species equation and calculates oxygen and hydrogen consumption rates. Figure 2 illustrates the geometry and computational cells of three flow patterns, including parallel, blocked bean-shaped, and unblocked bean-type. In the bean-shaped flow fields, pins are arranged in a consistent arrangement. This structure makes it possible to form parallel and spiral paths in the gas flow channels and achieve a more uniform flow distribution. In the bean-shaped flow field without blocking, the pins are placed, staggered, up to half of the channel height. However, in the blocked bean-shaped pattern, the height of pins and the channel is the same and pins are completely attached to the channel wall. This arrangement is employed for both cathode and anode sides. The geometrical features of flow patterns, along with physical and electrochemical properties, are presented in Table A2. other. The first reason is that the chemical and physical features of the materials vary and depend on the transfer parameters. The second reason is that conservation equations are coupled to each other using the Butler-Volmer equation. For instance, the Butler-Volmer equation is a function of hydrogen and oxygen concentrations; on the other hand, the species concentrations can be identified using the species equation. The Butler-Volmer equation enters as a source term into the species equation and calculates oxygen and hydrogen consumption rates. Figure 2 illustrates the geometry and computational cells of three flow patterns, including parallel, blocked bean-shaped, and unblocked bean-type. In the bean-shaped flow fields, pins are arranged in a consistent arrangement. This structure makes it possible to form parallel and spiral paths in the gas flow channels and achieve a more uniform flow distribution. In the bean-shaped flow field without blocking, the pins are placed, staggered, up to half of the channel height. However, in the blocked bean-shaped pattern, the height of pins and the channel is the same and pins are completely attached to the channel wall. This arrangement is employed for both cathode and anode sides. The geometrical features of flow patterns, along with physical and electrochemical properties, are presented in Table A2. Optimal meshes have a substantial impact on reducing computational cost and accuracy of results. Meshes with a small computational cell number decrease the solution accuracy. Theoretically, an increase in the number of computational cells makes the solution of governing equations more accurate; however, an excessive rise in the number of cells increases the computational cost, and enters the truncation error in the computational domain. The number of elements in the GDLs, CLs, and membrane is considered more than Optimal meshes have a substantial impact on reducing computational cost and accuracy of results. Meshes with a small computational cell number decrease the solution accuracy. Theoretically, an increase in the number of computational cells makes the solution of governing equations more accurate; however, an excessive rise in the number of cells increases the computational cost, and enters the truncation error in the computational domain. The number of elements in the GDLs, CLs, and membrane is considered more than other zones to improve computational accuracy because the electrochemical reactions occur in these regions. The numerical solution results perform after a grid-independent study on the PEMFC model. We used a personal computer with Intel(R) Core (TM) i7-6700 CPU @3.40 GHz 32 GB RAM to run the simulation. To make sure the results of the numerical modeling are independent of the computational grid size, about 593,920, 983,010, and 983,010 are considered for parallel, unblocked bean-shaped, blocked beantype, respectively. The number of elements in different zones of the PEMFC is presented in Table A3.

Results and Discussion
To verify the numerical results of the PEMFC model, a polarization curve of a singlecell PEMFC is compared to experimental measurements of Mazumder et al. [36], as exhibited in Figure 3a. Identical operating and geometrical conditions are performed for comparison, which includes the operating temperature of 50 • C, the membrane thickness of 0.23 mm, the working pressure of 1 atm, the GDL thickness of 0.254 mm, and the active area of 1.08 cm 2 . As is evident from this figure, the numerical results are comparable to the experimental measurements in these conditions, which are in good agreement. The maximum observed deviation between our modeling and experimental data is less than 4%. Hereupon, the numerical modeling can be utilized for a better understanding of complex electrochemical phenomena in the PEMFC. We also performed another comparison between our modeling and numerical study of Heidary et al. [15]. In this comparison, we simulated PEMFC with a pin flow field and reached a good agreement between the two models, as shown in Figure 3b. other zones to improve computational accuracy because the electrochemical reactions occur in these regions. The numerical solution results perform after a grid-independent study on the PEMFC model. We used a personal computer with Intel(R) Core (TM) i7-6700 CPU @3.40 GHz 32 GB RAM to run the simulation. To make sure the results of the numerical modeling are independent of the computational grid size, about 593,920, 983,010, and 983,010 are considered for parallel, unblocked bean-shaped, blocked beantype, respectively. The number of elements in different zones of the PEMFC is presented in Table A3.

Results and Discussion
To verify the numerical results of the PEMFC model, a polarization curve of a singlecell PEMFC is compared to experimental measurements of Mazumder et al. [36], as exhibited in Figure 3a. Identical operating and geometrical conditions are performed for comparison, which includes the operating temperature of 50 °C, the membrane thickness of 0.23 mm, the working pressure of 1 atm, the GDL thickness of 0.254 mm, and the active area of 1.08 cm 2 . As is evident from this figure, the numerical results are comparable to the experimental measurements in these conditions, which are in good agreement. The maximum observed deviation between our modeling and experimental data is less than 4%. Hereupon, the numerical modeling can be utilized for a better understanding of complex electrochemical phenomena in the PEMFC. We also performed another comparison between our modeling and numerical study of Heidary et al. [15]. In this comparison, we simulated PEMFC with a pin flow field and reached a good agreement between the two models, as shown in Figure 3b. It is essential to note that we select a voltage of 0.60 V as a corresponding point, which has a medium current density, and it is placed around the maximum power point of the PEMFC. Due to the higher diffusivity of hydrogen than air, the uniform reactant distribution at the cathode side is more critical than the anode side. Additionally, proper water management at the cathode side affects concentration drops and heat generation. Hence, the cathode side assessment is substantial for the appropriate design of the PEMFC, which will be discussed in the following. Figure 4 exhibits the streamlines at the cathode channel for the mentioned flow fields. It is found that the bean-shaped model without blockage has a better uniform reactant distribution than the other two models. In the blockage model, the height of pins is the same as channel height, which prevents passing the flow through pins. This fact forces the stream to cross around pins as a spiral flow. So, the vortices formation occurs behind the bean barriers. It is essential to note that we select a voltage of 0.60 V as a corresponding point, which has a medium current density, and it is placed around the maximum power point of the PEMFC. Due to the higher diffusivity of hydrogen than air, the uniform reactant distribution at the cathode side is more critical than the anode side. Additionally, proper water management at the cathode side affects concentration drops and heat generation. Hence, the cathode side assessment is substantial for the appropriate design of the PEMFC, which will be discussed in the following. Figure 4 exhibits the streamlines at the cathode channel for the mentioned flow fields. It is found that the bean-shaped model without blockage has a better uniform reactant distribution than the other two models. In the blockage model, the height of pins is the same as channel height, which prevents passing the flow through pins. This fact forces the stream to cross around pins as a spiral flow. So, the vortices formation occurs behind the bean barriers.   for parallel, unblocked bean-type, and blocked bean-shaped models, respectively. The particular arrangement of pins in the bean-type model without blockage causes a further increase in the flow velocity. In this flow field, the maximum velocity is observed in areas between two pins which are located adjacent to each other. These regions operate like a nozzle, which accelerates velocity as the path becomes narrower. However, the lowest velocity is identified behind pins. The flow has a high inertia force in these regions, and this can contribute to water aggregation. Subsequently, oxygen accessibility and the electrochemical reaction rate decrease in these areas. Moreover, gas mixture density reduces along the channel length due to oxygen consumption and conversion into water vapor. Since the density of water vapor is less than oxygen, a higher velocity is observed at the end of the channel at a constant mass flow rate.  for parallel, unblocked bean-type, and blocked bean-shaped models, respectively. The particular arrangement of pins in the bean-type model without blockage causes a further increase in the flow velocity. In this flow field, the maximum velocity is observed in areas between two pins which are located adjacent to each other. These regions operate like a nozzle, which accelerates velocity as the path becomes narrower. However, the lowest velocity is identified behind pins. The flow has a high inertia force in these regions, and this can contribute to water aggregation. Subsequently, oxygen accessibility and the electrochemical reaction rate decrease in these areas. Moreover, gas mixture density reduces along the channel length due to oxygen consumption and conversion into water vapor. Since the density of water vapor is less than oxygen, a higher velocity is observed at the end of the channel at a constant mass flow rate.  Figure 6 demonstrates the pressure contour of different flow fields along the cathode channel length for the stoichiometry of 20 and 30. The total pressure is the sum of two terms, i.e., dynamic pressure and static pressure. The static pressure term has a decreasing trend along the channel length because a part of the total energy of inlet gases is wasted to overcome losses, which arise from channel wall friction and the formation of small vortices behind the bean-shaped blockage. The gas mixture velocity along with the channel  Figure 6 demonstrates the pressure contour of different flow fields along the cathode channel length for the stoichiometry of 20 and 30. The total pressure is the sum of two terms, i.e., dynamic pressure and static pressure. The static pressure term has a decreasing trend along the channel length because a part of the total energy of inlet gases is wasted to overcome losses, which arise from channel wall friction and the formation of small vortices behind the bean-shaped blockage. The gas mixture velocity along with the channel length increases, and gas mixture density decreases. Based on the Bernoulli equation, the dynamic pressure term increases along the channel length; however, the portion of static pressure is greater than the dynamic pressure. Thus, it can be concluded that the total pressure decreases along the channel length. The highest amount of pressure drop is observed for the blocked pin-type at the cathode stoichiometric of 30. Because vortices created behind pins act as flow passing through pins, this increases the pressure drop in these areas and has a negative effect on fuel cell power output. The maximum amount of pressure drop at the cathode stoichiometric 30 for parallel, unblocking pin-type, and blocking pin-shaped models is 315.3 Pa, 11,844.4 Pa, and 20,762.7 Pa, respectively. This indicates that the unblocking pin-shaped pressure drop is nearly half of the blocking beanshaped flow field. Thus, this type of flow field has the potential to be used for PEMFC cooling with less pressure drop than the blocking-pin model at the stoichiometric of 30, but at a lower cathode stoichiometric (ξ cat = 20), these two models do not differ significantly in terms of pressure drop.  Figure 7 exhibits the distribution of oxygen concentration for the considered fl fields at the cathode GDL/CL interface. The comparison results show that bean-types h a higher oxygen concentration than the parallel flow channel due to more penetrat reactants into GDL. In bean-shaped flow channels, as air collides with pins, it devia toward the GDL and diffuses into it. As a result, the reactants' accessibility to the act surface area increases. The maximum oxygen diffusion into the cathode GDL at the ca ode stoichiometric of 30 in the unblocking pin-type, blocking bean-type, and parallel m els is 0.0065, 0.0069, and 0.0055, respectively. As oxygen availability within the CL creases, the electrochemical reaction rate and water production is raised, which leads improving fuel cell performance. Oxygen concentration is high at the entrance of the ca ode channel, and it gradually reduces along the channel length within electrochem reactions due to oxygen consumption. In areas below the channel's shoulder, less oxyg diffuses from the channel to the GDL due to more water accumulation. Since ribs are interface of BP and GDL, the accumulation of water impedes the transmission of electr current in these areas and reduces cell performance.  Figure 7 exhibits the distribution of oxygen concentration for the considered flow fields at the cathode GDL/CL interface. The comparison results show that bean-types have a higher oxygen concentration than the parallel flow channel due to more penetrating reactants into GDL. In bean-shaped flow channels, as air collides with pins, it deviates toward the GDL and diffuses into it. As a result, the reactants' accessibility to the active surface area increases. The maximum oxygen diffusion into the cathode GDL at the cathode stoichiometric of 30 in the unblocking pin-type, blocking bean-type, and parallel models is 0.0065, 0.0069, and 0.0055, respectively. As oxygen availability within the CL increases, the electrochemical reaction rate and water production is raised, which leads to improving fuel cell performance. Oxygen concentration is high at the entrance of the cathode channel, and it gradually reduces along the channel length within electrochemical reactions due to oxygen consumption. In areas below the channel's shoulder, less oxygen diffuses from the channel to the GDL due to more water accumulation. Since ribs are the interface of BP and surface area increases. The maximum oxygen diffusion into the cathode GDL at the ode stoichiometric of 30 in the unblocking pin-type, blocking bean-type, and parallel els is 0.0065, 0.0069, and 0.0055, respectively. As oxygen availability within the C creases, the electrochemical reaction rate and water production is raised, which lea improving fuel cell performance. Oxygen concentration is high at the entrance of the ode channel, and it gradually reduces along the channel length within electroche reactions due to oxygen consumption. In areas below the channel's shoulder, less ox diffuses from the channel to the GDL due to more water accumulation. Since ribs ar interface of BP and GDL, the accumulation of water impedes the transmission of elec current in these areas and reduces cell performance.   Figure 8 indicates water concentration distribution for the examined flow patterns at the cathode GDL/CL interface. The water concentration distribution has a similar trend to oxygen concentration; since oxygen is depleted along the channel length, water production decreases. The unblocking bean model has the highest amount of water production in comparison with other designs. In the stagnant areas of the PEMFC, both oxygen and water concentration are low, and the rate of electrochemical reactions decreases. In the parallel model, these areas can be observed at the side areas of the channel.  Figure 8 indicates water concentration distribution for the examined flow patterns at the cathode GDL/CL interface. The water concentration distribution has a similar trend to oxygen concentration; since oxygen is depleted along the channel length, water production decreases. The unblocking bean model has the highest amount of water production in comparison with other designs. In the stagnant areas of the PEMFC, both oxygen and water concentration are low, and the rate of electrochemical reactions decreases. In the parallel model, these areas can be observed at the side areas of the channel. A fuel cell needs a cooling system when it operates at a low stoichiometric ratio. For a small-sized fuel cell, the cooling target can only be achieved by raising the cathode stoichiometry ratio without requiring the cooling loop. In this regard, the cathode stoichiometry of 20 and 30 are utilized in this study. Figure 9 shows a comparison between the considered flow fields in terms of the temperature distribution in the middle of the cathode channel depth. As can be seen, the temperature along the flow channel rises due to heat transfer mechanisms, Joule heating, and the irreversibility heat. Lower temperature views are seen in bean-shaped models than the parallel design. This is because the convection heat transfer coefficient raises as velocity increases, and the heat transfer in the A fuel cell needs a cooling system when it operates at a low stoichiometric ratio. For a small-sized fuel cell, the cooling target can only be achieved by raising the cathode stoichiometry ratio without requiring the cooling loop. In this regard, the cathode stoichiometry of 20 and 30 are utilized in this study. Figure 9 shows a comparison between the Energies 2021, 14, 2494 14 of 23 considered flow fields in terms of the temperature distribution in the middle of the cathode channel depth. As can be seen, the temperature along the flow channel rises due to heat transfer mechanisms, Joule heating, and the irreversibility heat. Lower temperature views are seen in bean-shaped models than the parallel design. This is because the convection heat transfer coefficient raises as velocity increases, and the heat transfer in the pin vicinity improves. The highest temperature with the amount of 370 K is observed in the parallel flow field at the cathode stoichiometry of 20, which is above the desired temperature. According to the experimental work of Heras et al. [36], the optimum temperature of the fuel cell depends on the current, which can be calculated as below.
T opt = 0.53 I + 26.01 (40) Energies 2021, 14, x FOR PEER REVIEW 15 Based on Equation (40), the optimum temperature of PEMFC is 336 K at the cu of 70 A. Among the examined flow fields, the optimal performance of air-cooled is re to the unblocking bean-type at the stoichiometric of 20. The hottest point in unbloc bean-type reaches 337 K at the stoichiometry of 20, which is close to the desired temp ture of 336 K. It should be mentioned that all three models are almost similarly perform well in terms of thermal management at the cathode stoichiometry of 30.
The temperature distribution in the interface of GDL/CL on the cathode side is i trated in Figure 10. As depicted in the figure, using bean-type models leads to mor vorable cooling than the parallel design. Pins, which are half of the channel height i unblocked bean-shape, play the role of the fin. In other words, pins that are up to middle of the channel height and are not in contact with the GDL have a cooling ro air flows through them. This causes the conduction heat transfer to increase; conseque the temperature further drops in these regions. Further, increasing velocity in pinmodels improves the PEMFC cooling due to the improvement of the heat conve mechanism. Based on Equation (40), the optimum temperature of PEMFC is 336 K at the current of 70 A. Among the examined flow fields, the optimal performance of air-cooled is related to the unblocking bean-type at the stoichiometric of 20. The hottest point in unblocking bean-type reaches 337 K at the stoichiometry of 20, which is close to the desired temperature of 336 K. It should be mentioned that all three models are almost similarly performing well in terms of thermal management at the cathode stoichiometry of 30.
The temperature distribution in the interface of GDL/CL on the cathode side is illustrated in Figure 10. As depicted in the figure, using bean-type models leads to more favorable cooling than the parallel design. Pins, which are half of the channel height in the unblocked bean-shape, play the role of the fin. In other words, pins that are up to the middle of the channel height and are not in contact with the GDL have a cooling role as air flows through them. This causes the conduction heat transfer to increase; consequently, the temperature further drops in these regions. Further, increasing velocity in pin-type models improves the PEMFC cooling due to the improvement of the heat convection mechanism. trated in Figure 10. As depicted in the figure, using bean-type models leads to more favorable cooling than the parallel design. Pins, which are half of the channel height in the unblocked bean-shape, play the role of the fin. In other words, pins that are up to the middle of the channel height and are not in contact with the GDL have a cooling role as air flows through them. This causes the conduction heat transfer to increase; consequently, the temperature further drops in these regions. Further, increasing velocity in pin-type models improves the PEMFC cooling due to the improvement of the heat convection mechanism.   Figure 11 displays the current density distribution at the interface of the cathode GDL/CL. Under rib areas in the parallel flow field and regions under bean-pins, the current density is lower than other regions because of the water accumulation under ribs and contact resistance between BP and GDL interface. Ribs are the interface of BP and GDL, and the electrical current flows in these places. The more accumulated water under rib areas, the harder transferring the electrical current becomes. A high amount of current density is observed for the margins of the cell due to the lake of the accumulated water in these areas. It can be inferred that the distribution of current density in the parallel model is more influenced by concentration losses, which arise from insufficient access of reactant to the reaction surface. The current density decreases in the middle of the unblocking bean-type due to the lack of direct contact of pins with the GDL. As depicted in Figure 11, the maximum current density of blockage and non-blockage pin-design at the cathode stoichiometry of 20 are 5.9 A cm −2 and 4 A cm −2 , respectively. Although the rate of electrochemical reaction and the extraction of the current density increases at a higher stoichiometry ratio, increasing the cathode stoichiometry ratio of more than 20 doses not affect the electrochemical reaction rate. A higher cathode stoichiometric ratio than 20 can only facilitate fuel cell cooling, as illustrated in Figure 9.  Figure 11 displays the current density distribution at the interface of the cathode GDL/CL. Under rib areas in the parallel flow field and regions under bean-pins, the current density is lower than other regions because of the water accumulation under ribs and contact resistance between BP and GDL interface. Ribs are the interface of BP and GDL, and the electrical current flows in these places. The more accumulated water under rib areas, the harder transferring the electrical current becomes. A high amount of current density is observed for the margins of the cell due to the lake of the accumulated water in these areas. It can be inferred that the distribution of current density in the parallel model is more influenced by concentration losses, which arise from insufficient access of reactant to the reaction surface. The current density decreases in the middle of the unblocking bean-type due to the lack of direct contact of pins with the GDL. As depicted in Figure 11, the maximum current density of blockage and non-blockage pin-design at the cathode stoichiometry of 20 are 5.9 A cm −2 and 4 A cm −2 , respectively. Although the rate of electrochemical reaction and the extraction of the current density increases at a higher stoichiometry ratio, increasing the cathode stoichiometry ratio of more than 20 doses not affect the electrochemical reaction rate. A higher cathode stoichiometric ratio than 20 can only facilitate fuel cell cooling, as illustrated in Figure 9.  Figure 12 illustrates the polarization curves of three considered flow fields at the cathode stoichiometric ratios of 20 and 30. As can be seen, the worst performance is attained in the blocking pin-shaped model due to lower output power in comparison to other flow fields. There is not much difference between the bean-shaped design without  Figure 12 illustrates the polarization curves of three considered flow fields at the cathode stoichiometric ratios of 20 and 30. As can be seen, the worst performance is attained in the blocking pin-shaped model due to lower output power in comparison to other flow fields. There is not much difference between the bean-shaped design without blocking and parallel model under medium current densities. The unblocked bean-type design shows a better performance than parallel channels at high current densities where the concentration losses are the main reason for the voltage drop. This better performance can be attributed to the better mass transfer phenomena and the distribution of the proper reactants on the catalyst surface. The influence of oxygen consumption is more pronounced at high current densities because oxygen demand is high in this range; hence, the curves for considered models are significantly different at high current densities. Oxygen depletion before reaching to entire catalyst surface is the reason for limiting the current density at the blocking bean-shaped. The current density at the voltage of 0.6 V and cathode stoichiometric of 20 is 0.913 A cm −2 , 0.974 A cm −2 , and 0.683 A cm −2 , for parallel, unblocking, and blocking pin-type models, respectively. The fuel cell performance improves at a higher stoichiometry because of more oxygen access to the reaction surface, which results in raising the electrochemical reaction rate and generating more current. However, at a higher stoichiometric ratio (more than 20), there is no noticeable difference in system performance. The current density of parallel, unblocking, and blocking pin-type patterns is 0.940 A cm −2 , 0.998 A cm −2 , and 0.693 A cm −2 , respectively, at the voltage of 0.6 V and cathode stoichiometric of 30. which results in raising the electrochemical reaction rate and generating more current. However, at a higher stoichiometric ratio (more than 20), there is no noticeable difference in system performance. The current density of parallel, unblocking, and blocking pin-type patterns is 0.940 A cm −2 , 0.998 A cm −2 , and 0.693 A cm −2 , respectively, at the voltage of 0.6 V and cathode stoichiometric of 30.  Figure 13 indicates the power density curve for three flow fields. The highest power density for the blocking, unblocking bean-types, and parallel design are 0.564 W cm −2 , 0.599 W cm −2 , and 0.416 W cm −2 at the operating voltage of 0.6 and cathode stoichiometric of 30. By increasing the cathode stoichiometry from 20 to 30, the maximum power density increases by almost 2.91%, 2.56%, and 1.46% for parallel, unblocking pin-type, and blocking pin-shaped models, respectively. 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 1.3 1 Figure 13 indicates the power density curve for three flow fields. The highest power density for the blocking, unblocking bean-types, and parallel design are 0.564 W cm −2 , 0.599 W cm −2 , and 0.416 W cm −2 at the operating voltage of 0.6 and cathode stoichiometric of 30. By increasing the cathode stoichiometry from 20 to 30, the maximum power density increases by almost 2.91%, 2.56%, and 1.46% for parallel, unblocking pin-type, and blocking pin-shaped models, respectively. Previous research has shown that less than 3 kW power is needed for UAV applications, which can be supplied using a PEMFC stack [37]. Recently, UK-based Intelligent Energy launched a UAV powered by a 2.4 kW fuel cell stack. They mentioned that more flight time (2 h, i.e., more than five times) was achieved using fuel cells over batteries (25 min). The company utilized a 2.4 kW fuel cell stack for the UAV with 11 L hydrogen, which was stored in the high-pressure hydrogen tank [37]. In the present study, the target is providing the power of 2.5 kW for a UAV using the PEMFC stack. Based on the presented numerical results, the mass and volume characteristics of a single fuel cell for different flow fields are shown in Table 1. Asghari et al. [38] reported that the optimum endplate thickness of the PEMFC stack with an output power of 5 kW is 35 mm. The mass and volume of endplates were achieved at 1.89 kg and 0.7 L, respectively, based on the active area of 100 cm 2 and using aluminum as the endplate's material. We use the properties of mentioned endplates in the fuel cell stack design. Given that the active area in numerical modeling is considered small due to computational constraints and is different from their examined active area, we get larger the scale of the PEMFC stack for the design of the PEMFC stack. Accordingly, the characteristics of the PEMFC stack, including the number of cells, volume, mass, and parasitic

Current density[A cm -2 ]
Power density[W cm Previous research has shown that less than 3 kW power is needed for UAV applications, which can be supplied using a PEMFC stack [37]. Recently, UK-based Intelligent Energy launched a UAV powered by a 2.4 kW fuel cell stack. They mentioned that more flight time (2 h, i.e., more than five times) was achieved using fuel cells over batteries (25 min). The company utilized a 2.4 kW fuel cell stack for the UAV with 11 L hydrogen, which was stored in the high-pressure hydrogen tank [37]. In the present study, the target is providing the power of 2.5 kW for a UAV using the PEMFC stack. Based on the presented numerical results, the mass and volume characteristics of a single fuel cell for different flow fields are shown in Table 1. Asghari et al. [38] reported that the optimum endplate thickness of the PEMFC stack with an output power of 5 kW is 35 mm. The mass and volume of endplates were achieved at 1.89 kg and 0.7 L, respectively, based on the active area of 100 cm 2 and using aluminum as the endplate's material. We use the properties of mentioned endplates in the fuel cell stack design. Given that the active area in numerical modeling is considered small due to computational constraints and is different from their examined active area, we get larger the scale of the PEMFC stack for the design of the PEMFC stack. Accordingly, the characteristics of the PEMFC stack, including the number of cells, volume, mass, and parasitic power consumption due to pressure drop in cathode channels, are stated in Table 2. It is worth noting that the maximum output power is obtained around an operating voltage of 0.60 V.
The unblocking bean-type generates more power than other models and needs fewer cells number for producing an output power of 2.5 kW. The number of required cells to generate the net power of 2.5 kW for the parallel, unblocking, and blocking pin-shaped models is 60, 55, and 80, respectively. Therefore, unblocking bean-shaped has a significant influence on the reduction of PEMFC weight and volume. The power consumption to overcome the pressure drop along the cathode gas channel is calculated by the following equation: The power consumption of regenerative compressor for providing the required air and dominating pressure drop within cathode channels is 0.84 W, 29.8 W, and 69.9 W for parallel, unblocking bean-shaped, blocking bean-shaped, respectively. The PEMFC stack operating voltages for parallel, unblocking, and blocking pin models are 36 V, 33 V, and 48 V, respectively. Based on the operating voltage (33-48 V) of the fuel cell stack for different flow fields and pressure drop along the cathode channels, we select a spark-free regenerative compressor (VRB8-25) from VAIREX company with a proper pressure ratio to supply the required air [39]. This regenerative compressor has a weight of 2 kg, and its operating voltage is in the range of 24-48 V, which is within the operating voltage range of the fuel cell, and we do not need to use a converter [39]. Figure 14 represents a comparison between three flow fields in terms of mass and volume power density. The maximum amount of volume power density, i.e., 1.1 kW L −1 , is related to the bean-shaped without blocking, and the minimum amount of power density is allocated to bean-type with blocking. In terms of mass power density, the maximum value of 0.2 kW kg −1 corresponds to the unblocked bean flow field. power consumption due to pressure drop in cathode channels, are stated in Table 2. It is worth noting that the maximum output power is obtained around an operating voltage of 0.60 V. The unblocking bean-type generates more power than other models and needs fewer cells number for producing an output power of 2.5 kW. The number of required cells to generate the net power of 2.5 kW for the parallel, unblocking, and blocking pinshaped models is 60, 55, and 80, respectively. Therefore, unblocking bean-shaped has a significant influence on the reduction of PEMFC weight and volume. The power consumption to overcome the pressure drop along the cathode gas channel is calculated by the following equation: The power consumption of regenerative compressor for providing the required air and dominating pressure drop within cathode channels is 0.84 W, 29.8 W, and 69.9 W for parallel, unblocking bean-shaped, blocking bean-shaped, respectively. The PEMFC stack operating voltages for parallel, unblocking, and blocking pin models are 36 V, 33 V, and 48 V, respectively. Based on the operating voltage (33-48 V) of the fuel cell stack for different flow fields and pressure drop along the cathode channels, we select a spark-free regenerative compressor (VRB8-25) from VAIREX company with a proper pressure ratio to supply the required air [39]. This regenerative compressor has a weight of 2 kg, and its operating voltage is in the range of 24-48 V, which is within the operating voltage range of the fuel cell, and we do not need to use a converter [39]. Figure 14 represents a comparison between three flow fields in terms of mass and volume power density. The maximum amount of volume power density, i.e., 1.1 kW L −1 , is related to the bean-shaped without blocking, and the minimum amount of power density is allocated to bean-type with blocking. In terms of mass power density, the maximum value of 0.2 kW kg −1 corresponds to the unblocked bean flow field. In summary, it can be concluded that the PEMFC stack with the unblocked bean-type has the best performance in terms of temperature distribution and mass power density. In summary, it can be concluded that the PEMFC stack with the unblocked bean-type has the best performance in terms of temperature distribution and mass power density. This novel flow pattern is suggested for aerial applications that aim to reduce weight and attain optimal heat management.

Conclusions
Energy density is the critical factor affecting the endurance of UAVs' flight time. Currently, lithium-ion batteries are used as an energy source for UAVs and have a small energy density, which causes short flight times. It is necessary to increase the flight time of UAVs using alternative power sources such as fuel cells. In recent years, PEMFCs have been considered, and they have the potential to be used as a power source of UAVs, but there should be improvements in their power density. The weight and volume of the fuel cell stack can be significantly reduced by improving the configuration and the embedded flow field within BPs. As such, the proper design of PEMFC flow fields contributes substantially to increasing power density and UAV flight times. In this regard, a comprehensive study is carried out on the three flow fields, i.e., parallel, bean-shaped without blocking, and beanshaped with blocking, for aerial applications in terms of mass and volume power density. Another goal of this study was to examine the cooling performance of the considered flow fields to maintain optimum operating temperature during PEMFC operation. The key findings of this paper are summarized below: 1.
The specific arrangement of pins in the unblocking bean-type allows flow velocity to increase further. The highest velocity is found in this flow field in areas between two pins, where they act as a nozzle and increase velocity as the path becomes narrow.

2.
For blockage pin-type, the highest amount of pressure drop is observed for a cathode stoichiometric of 30 and is nearly 20,762.7 Pa due to creating vortices behind pins. However, there is not a notable difference between bean-type models in terms of pressure drop in a lower cathode stoichiometric (ξ cat = 20).

3.
In pin-type models, oxygen penetration into the GDL is more than a parallel design due to the specific arrangement of pins. When air collides with beans, it diverts and diffuses to the GDL. This increases the accessibility of the reactants to the active surface area.

4.
In bean-shaped models, a slight difference in temperature is seen in the space between two pins because the coefficient of convection heat transfer raises as velocity increases in these areas, which can facilitate heat management.

5.
Pins that are placed until the middle of the channel's depth play the role of a fin in the unblocked bean-shape. This leads to enhancing heat conduction and improving the PEMFC cooling using a further drop in temperature. The hottest point in the bean-shape without blocking in the stoichiometry of 20 is 337.9 K, which is close to the optimal temperature of 336 K. With increasing stoichiometry from 20 to 30, the maximum temperature for parallel, unblocking, and blocking pin-type models reduce from 370, 337.9, and 341.6 to 344.6, 331.5, and 332.3, respectively. By using an unblocked bean-shape in a stoichiometry of 20, we achieve optimal cooling and do not need to increase cathode stoichiometry for the cooling target. For other flow fields, we need to increase the cathode stoichiometry to reach the desired temperature. 6.
The current density distribution becomes more uniform as the cathode stoichiometry increases. That causes the electrochemical reaction rate and the electrical current density to increase. As the reactant diffuses into active surface areas, and the electrochemical reaction occurs, the current density reduces along the active surface area due to oxygen consumption. Nonetheless, there is no considerable difference in system performance at a higher stoichiometric ratio (more than 20). 7.
The highest amount of power density is achieved in the bean-type without blockage at a specific voltage of 0.60, which is about 6.2% and 43.9% higher than the parallel and blocking bean-type, respectively. Therefore, the unblocked bean-shape produces more power for a given active area and requires fewer cells to generate power. For the parallel, blocking, and unblocking pin models, the number of needed cells to generate net power of 2.5 kW is 60, 80, and 55, respectively. 8.
Among different investigated flow fields, the bean-shape without blockage is an ideal flow field in the PEMFC for aerial applications with the maximum amount of volume power density, i.e.,1.1 kW L −1 , and maximum mass power density of 0.2 kW kg −1 .

Conflicts of Interest:
The authors declare no conflict of interest.
Appendix A Table A1. Considered boundary condition in the numerical modeling of PEMFC.