PRECISION DESIGN OF PROTON EXCHANGE MEMBRANE FUEL CELL FLOW FIELD BASED ON WATER AND HEAT TRANSFER MECHANISM

Abstract

The article establishes a three-dimensional multiphase proton exchange membrane single-cell model and investigates the impact of precision flow field design on the electrochemical characteristics, heat mass transfer properties, and phase change characteristics of fuel cells. The simulation model is analyzed using COMSOL 6.0 multiphysics software and validated using experimental data under the same operating conditions. The research results indicate that precision flow field design can enhance the electrochemical characteristics of proton exchange membrane fuel cells (PEMFC) by reducing concentration overpotentials through improved gas mass transfer, thus increasing the performance of cell. Precision flow field design improves gas and current density distribution uniformity, promoting more uniform electrochemical reactions and enhancing cell durability. Finally, precision flow field design increases gas mass transfer rates, effectively removing liquid water from the interior and facilitating smoother gas transport to the active regions. This study provides new insights into flow field design for fuel cells.

Keywords:
fuel cell; flow field; precision design; water heat transfer; phase change

INTRODUCTION

Proton exchange membrane fuel cells (PEMFCs), with their high energy density, no air pollution, and high energy conversion efficiency, have emerged as one of the ideal power sources in the 21^st century. However, the transient response characteristics of fuel cells remain a critical factor hindering their commercialization. Recognizing this, the academic community has employed various methods to enhance their performance, and flow field optimization is one of these methods.¹, ², ³ The flow field on the bipolar plates plays a crucial role in PEMFC operation. It evenly distributes the fuel and oxidant required for electrochemical reactions to ensure uniform current density distribution. It facilitates the smooth removal of reaction-generated water under the influence of exhaust gases and entrainment.⁴ The focus of flow field optimization design primarily revolves around flow field geometry and dimensions.

Traditional flow field design: currently, commonly used flow field structures include traditional ones like parallel flow fields and serpentine flow fields, as well as newer designs such as 3D flow fields, biomimetic flow fields, and radial flow fields. In the case of parallel flow fields, gas experiences low flow resistance when flowing through it. However, these flow fields need better water management characteristics. Water vapor generated during electrochemical reactions tends to condense within the channels, blocking the gas transport and resulting in uneven distribution of reaction gases in the electrocatalytic region. To address these issues, researchers have made various attempts.⁵, ⁶, ⁷ Timurkutluk and Chowdhury⁸ investigated the impact of variable cross-sectional structures with width and height gradually narrowing along the length on straight flow fields. The results confirmed that such variable cross-sectional structures can effectively improve the distribution of oxygen and water on the cathode catalyst layer.

Serpentine flow field structure is a structure widely used in recent years. It is a continuous channel of reciprocating and bending from the inlet to the outlet. Compared with the parallel flow field, the serpentine flow field has better water removal performance, and it is not easy to appear flooding in the flow field. However, because the serpentine flow field has many inflection points that can change the gas flow direction sharply, the flow field easily produces a large pressure drop, resulting in an insufficient reaction gas supply in the second half of the flow field. Some experts⁹, ¹⁰ proposed multiple snake flow fields to solve the problem of excessive pressure drop in the snake flow field. The simulation analysis shows that the PEMFC performance can be improved by adjusting multiple serpentine flow fields’ number, size, and length.

New flow field design: based on the traditional flow field, some researchers proposed a new flow field, such as a 3D flow field, bionic flow field, and radial flow field, and conducted simulation research.

Toyota proposed the 3D flow field and has now achieved mass production on its hydrogen cell model Mirai. The structural characteristic of a 3D flow field is that many microfow fields are processed based on the structure of the straight flow field on the cathode side. The air spreads through the microcurrent field and can flow more fully to the catalyst layer to participate in the reaction. At the same time, the product water can also be quickly removed from the flow field through the microcurrent field to avoid the problem of water aggregation and blocking the flow field. Niu et al.¹¹ designed two new 3D flow fields of cathode flow field by increasing the baffle at the surface to improve the gas diffusion layer and enhance the mass transfer between the flow field and the gas diffusion layer. Zhang et al.¹² established a three-dimensional multiphase numerical model of PEMFC and studied PEMFC with a 3D fine pore flow field. They found that a 3D thin pore flow field can significantly improve the reactive gas supply from the flow field to the porous electrode and simultaneously promote the removal of liquid water in PEMFC.

The biomimetic flow field draws on biological structures such as leaf veins and animal lungs, and the advantages of such structures in mass transfer and heat transfer have been verified in relevant literature.¹³, ¹⁴ Kang et al.¹⁵ studied the gas transmission phenomenon of the wing vein flow field. Compared with the serpentine flow field and the leaf vein flow field, the wing vein flow field can significantly improve the uniformity of the reaction gas concentration distribution in the gas diffusion layer (GDL), which is more conducive to strengthening the mass transfer of the reaction gas. Damian-Ascencio et al.¹⁶ studied the role of mass transfer reinforcement of tree-like flow field structure in the design of the introduced cross-flow field. The results show that introducing a tree-type flow field can better reduce the water content in the proton exchange membrane and improve water transmission performance.

The radial flow field structure adopts the central intake arrangement, and the gas flows from the inlet at the center of the flow field to the periphery of the flow field and finally exits the PEMFC through the outlet distributed at the edge of the circular bipolar plate. Compared with the two new flow fields above, the radial flow field has the advantages of a short intake path, easy processing, and convenient stacking under the same working area. Cano-Andrade¹⁷ selected a quarter of the radial flow field region and simulated the current density distribution characteristics of PEMFC and the mass transfer characteristics of GDL with the structure of 4, 8, and 12 flow fields. The results show that the three factors significantly improve performance than the traditional parallel and snake flow fields.

Optimized design of the flow field size: after many experiments and simulation studies, many researchers¹⁸, ¹⁹ believe that the length of the flow field greatly impacts cell performance. These effects are mainly manifested in the fact that the longer the flow field, the greater the pressure loss, the lower the gas concentration of the reaction in the later stage, and the flooding phenomenon, thus reducing the cell performance. There is little research on the impact of flow field depth on cell performance, but many designers believe that a shallow flow field can get a signifcant flow rate and prevent flooding to get a better performance.²⁰, ²¹ There are many studies²², ²³ on the influence of flow field width on PEMFC performance, and many experimental and simulation results show that a small ridge width and a large flow field width can promote the reaction gas mass transfer, thus improving cell performance. Most researchers²⁴, ²⁵, ²⁶, ²⁷, ²⁸ are studying the section shape problem from the perspective of pressure drop loss and processing methods. For example, the triangular and semicircular section is more conducive to reducing the pressure drop than the rectangular section; the increase of the slope angle will cause the increase of the maximum power density and current density of PEMFC, and the trapezoidal section flow field can effectively remove liquid water and so on. Some scholars²⁹ believe that changing the internal structure of the flow field, such as adding a baffle inside the flow field, can effectively improve the quality transfer performance of PEMFC.

In conclusion, the following problems exist in the flow field optimization design: (i) there are many studies on the flow field form and the optimization design of the flow field structure size. The new flow field design is constantly proposed, and the relationship between the flow field size and PEMFC performance is relatively clear. However, there are few studies on more precise, small-size flow field design; (ii) to study the optimal design of the flow field, the evaluation indexes are often external characterization, such as current density and power density, as well as internal characteristics, such as gas and electric density distribution. The influence of the change of the flow field on the mass transfer and water vapor phase transition inside the cell is ignored.

Given this, this paper proposes a flow field precision design strategy to compare and analyze the infuence of flow field forms with different size accuracy, such as width × height: 1 mm × 1 mm, 0.7 mm × 0.7mm, and 0.5 mm × 0.5 mm, on PEMFC performance. The influence of the precision design of the exploration field on the water heat transfer and water vapor phase transition law of the fuel cell provides a new idea for the flow field design of the fuel cell flow field.

EXPERIMENTAL

Methodology

Mathematical equations

Fuel cell internal includes electrochemical reaction, phase change, mass transfer, heat transfer, and other complex processes. Therefore, in the simulation calculation, it is necessary to use multiple mathematical equations to make the calculation more accurate. The main mathematical equation contains the following categories.

Electrochemical equations

The electrochemical equations reveal the crucial electrochemical reactions occurring in fuel cells. Typical fuel cell reactions involve the electrochemical oxidation-reduction processes of hydrogen and oxygen at the anode and cathode, respectively. These equations provide the foundation for understanding the energy conversion mechanisms within fuel cells. The equations are as follows:

(i) Exchange current density

(1) $i_{0,a}=i_{0,a}^{ref}\left(1-s_{ice}-s_l\right){\left(\frac{c_H{}_{{}_2}}{c_{H_2}^{ref}}\right)}^{0.5}\times \exp \left(-1400\left(\frac{1}{T}-\frac{1}{298.15}\right)\right)$

where i_0a is the anodic exchange current density; $i_{0,a}^{ref}$ is the anode reference exchange current density; s_l is the saturation of liquid water; s_ice is the saturation of ice; $c_{H_2}$ is the hydrogen concentration; $c_{H_2}^{ref}$ is the hydrogen reference concentration; T is the cell temperature.

(2) $i_{0,a}=i_{0,a}^{ref}\left(1-s_{ice}-s_l\right){\left(\frac{c_{H_2}}{c_{H_2}^{ref}}\right)}^{0.5}\times \exp \left(-7900\left(\frac{1}{T}-\frac{1}{298.15}\right)\right)$

(ii) Open circuit voltage

(3) $V_{\text{out}}=V_{nernst}+{\text{η}}_{act,a}+{\text{η}}_{act,c}+{\text{η}}_{conc,a}+{\text{η}}_{conc,c}+{\text{η}}_{ohmic}$

where V_out is the output voltage; V_nernst is the Nister voltage; η_act,a is the anodic activation of the overpotential; η_act,c is the cathodic activation by overpotential; η_conc,a is the anode concentration difference and overpotential; η_conc,c is the cathodic concentration overpotential; η_ohmic is the ohm overpotential.

(iii) Nernst voltage

(4) $V_{nernst}=1.229-0.846\times {10}^{-3}\times \left(T-298.15\right)+\frac{RT}{2F}\ln \left(p_{H_2}{p}_{O_2}^{0.5}\right)$

where p_H2 is the hydrogen pressure; p_O2 is the oxygen partial pressure.

(iv) Activated over-potential

(5) ${\text{η}}_{act,a}=-\frac{RT}{2{\alpha }_{a}F}\ln \left(\frac{I}{{\delta }_{CL,a}{i}_{0,a}}\right)$

(6) ${\text{η}}_{act,c}=-\frac{RT}{4{\alpha }_{c}F}\ln \left(\frac{I}{{\delta }_{CL,c}{i}_{0,c}}\right)$

where α_a and α_c is the transfer coefficient of the anode and the cathode; I is current density; δ_CL is catalytic layer (CL) thickness.

(v) Concentration difference over potential

(7) ${\text{η}}_{conc,a}=-\frac{RT}{2{\alpha }_{a}F}\ln \left(1-\frac{I}{I_{D,a}}\right)$

(8) ${\text{η}}_{conc,c}=-\frac{RT}{4{\alpha }_{a}F}\ln \left(1-\frac{I}{I_{D,c}}\right)$

where I_D,a and I_D,c are the limiting current density of the anode and the cathode, respectively.

(vi) Limiting current density

(9) $I_{D,a}=\frac{2F{c}_{H_2}}{\frac{{\delta }_{CL,a}}{2D_{H_2,CL}^{eff}}+\frac{{\delta }_{GDL}}{D_{H_2,GDL}^{eff}}}$

(10) $I_{D,c}=\frac{4F{c}_{O_2}}{\frac{{\delta }_{CL,c}}{2D_{O_2,CL}^{eff}}+\frac{{\delta }_{GDL}}{D_{O_2,GDL}^{eff}}}$

(11) $D_i^{eff}=D_i\left(1-s_l\right){\epsilon }^{1.5}$

where δ_GDL is GDL thickness; D^eff is the reference diffusion coefficient; ε is porosity.

(vii) Coefficient of diffusion

(12) $D_{H_2}=1.055\times {10}^{-4}{\left(\frac{T}{333.15}\right)}^{1.5}$

(13) $D_{O_2}=2.652\times {10}^{-5}{\left(\frac{T}{333.15}\right)}^{1.5}$

(14) $D_{vp}=2.982\times {10}^{-5}{\left(\frac{T}{333.15}\right)}^{1.5}$

(15) $D_{N_2}=2.652\times {10}^{-5}{\left(\frac{T}{333.15}\right)}^{1.5}$

(viii) Ohm overpotential

(16) ${\text{η}}_{ohmic}=-\left(\frac{{\delta }_{GDL}}{{\sigma }_{s,GDL}}+\frac{{\delta }_{CL}}{{\sigma }_{s,CL}}+\frac{{\delta }_{MEM}}{{\sigma }_{MEM}}\right)$

where δ_MEM is the thickness of the membrane; σ_s,GDL, σ_s,CL, σ_s,MEM are the electrical conductivity of the GDL, CL and membrane.

Conservation equations

In fuel cells, the conservation equations for energy, mass, and charge describe the transfer and conversion of these quantities within the system. These equations are crucial for understanding the operational principles of fuel cells and for optimizing their design and performance. Here are some important conservation equations in fuel cells and their roles:

(i) Mass conservation

(17) $\frac{\partial }{\partial t}\left(\epsilon \left(1-s_l\right){\rho }_g\right)+\nabla \times \left({\rho }_g{\overrightarrow{u}}_g\right)=S_m$

(18) ${\rho }_g=p_g{\left(RT{\displaystyle \sum _i\frac{Y_i}{M_i}}\right)}^{-1}$

where ρ_g is the gas mixing density; Swhere ρ_g is the gas mixing density; S_m is the quality source phase; Y is the mole fraction; M is the relative molecular weight.

where ρ_g is the gas mixing density; S_m is the quality source phase; Y is the mole fraction; M is the relative molecular weight.

_m is the quality source phase; Y is the mole fraction; M is the relative molecular weight.

(ii) Momentum conservation

(19) $\frac{\partial }{\partial t}\left(\frac{{\rho }_g{\overrightarrow{u}}_g}{\epsilon \left(1-s_l\right)}\right)+\nabla \left(\frac{{\rho }_g{\overrightarrow{u}}_g{\overrightarrow{u}}_g}{{\epsilon }^2{\left(1-s_l\right)}^2}\right)=-\nabla p_g+{\mu }_g\nabla \left(\nabla \left(\frac{{\overrightarrow{u}}_g}{\epsilon \left(1-s_l\right)}\right)+\nabla \left(\frac{{\overrightarrow{u}}_g^T}{\epsilon \left(1-s_l\right)}\right)\right)-\frac{2}{3}{\mu }_g\nabla \left(\nabla \left(\frac{{\overrightarrow{u}}_g}{\epsilon \left(1-s_l\right)}\right)\right)+S_u$

where S_u is the momentum source phase.

(iii) Material conservation (hydrogen gas, oxygen, and water vapor)

(20) $\frac{\partial }{\partial t}\left(\epsilon \left(1-s_l\right){\rho }_g\right)+\nabla \left({\rho }_g{\overrightarrow{u}}_g{Y}_i\right)=\nabla \left({\rho }_g{D}_i^{eff}\nabla Y_i\right)+S_i$

(21) ${\mu }_g={\displaystyle \sum _i\frac{X_i{\mu }_i}{{\displaystyle {\sum }_j{X}_i{\psi }_{\mu ij}}}}$

(22) ${\psi }_{\mu ij}=\frac{{\left[1+{\left(\frac{{\mu }_i}{{\mu }_j}\right)}^{0.5}{\left(\frac{M_i}{M_j}\right)}^{0.25}\right]}^{0.5}}{{\left[8\left(1+\frac{M_i}{M_j}\right)\right]}^{0.5}}$

where µ_g is the ideal mixed gas dynamic viscosity; X is the mole fraction; µ_i is the dynamic viscosity of the gas; ψ is the relative molar mass; $D_i^{eff}$ is the gas reference diffusivity factor; S_i is the material source phase.

(iv) Gas power viscosity

(23) ${\mu }_{H_2}=3.205\times {10}^{-3}{\left(\frac{T}{293.85}\right)}^{1.5}{\left(T+72\right)}^{-1}$

(24) ${\mu }_{O_2}=8.46\times {10}^{-3}{\left(\frac{T}{292.25}\right)}^{1.5}{\left(T+127\right)}^{-1}$

(25) ${\mu }_{vp}=7.512\times {10}^{-3}{\left(\frac{T}{291.15}\right)}^{1.5}{\left(T+120\right)}^{-1}$

(v) Liquid water conservation

(26) $\frac{\partial \left(\epsilon s_l{\rho }_l\right)}{\partial t}+\nabla \left({\rho }_l{D}_l\nabla s_l\right)+\nabla \iota {\rho }_l{\overrightarrow{u}}_g=S_l$

(27) $\iota =\frac{K_l{\mu }_g}{K_g{\mu }_l}\text{iota}$

(28) ${\mu }_l=2.414\times {10}^{-5}\times {10}^{\frac{247.8}{T-140}}$

(29) $D_l=-\frac{K_l}{{\mu }_l}\frac{d{p}_c}{d{s}_l}$

(30) $K_g=K_0{\left(1-s_l\right)}^3$

(31) $K_l=K_0{s}_l^3$

(32) $p_c=F_{st}\cos \theta \frac{\epsilon }{K_o}J\left(s\right)$

(33) $J\left(s\right)=\{\begin{array}{cc}1.417\left(1-s_l\right)-2.120{\left(1-s_l\right)}^2+1.263{\left(1-s_l\right)}^3& ,\theta <90°\\ 1.417s_l-2.120s_l^2+1.263s_l^3& \theta \ge 90°\end{array}$

(34) $\begin{array}{cc}{F}_{st}=-0.0001676T+0.1218& 273.15K\le T\le 373.15K\end{array}$

where ρ_l is the density of water; µ_l> is the liquid hydrodynamic viscosity; ι is the interface resistance coefficient; D_l is the liquid-water diffusion coefficient; p_c is the capillary pressure; F_st is the surface tension between liquid water and gas; K₀ is the permeability; θ is the contact angle; K_l is the effective liquid permeability; K_g is the relative gas-phase permeability.

(vi) Ionomer hydration water conservation

(35) $\frac{{\rho }_{mem}}{EW}\frac{\partial \left(\omega {\lambda }_{nf}\right)}{\partial t}-\frac{{\rho }_{mem}}{EW}\nabla \left({\omega }^{1.5}{D}_{nmw}\nabla {\lambda }_{nf}\right)=S_{nmw}$

(36) $D_{nmw}=\{\begin{array}{cc}3.1\times {10}^{-7}{\lambda }_{nf}\left(\exp \left(0.28{\lambda }_{nf}\right)-1\right)\exp \left(\frac{-2346}{T}\right)& ,0<{\lambda }_{nf}<3\\ 4.17\times {10}^{-8}{\lambda }_{nf}\left(161\exp \left(-{\lambda }_{nf}\right)+1\right)\exp \left(\frac{-2346}{T}\right)& 3\le {\lambda }_{nf}<17\end{array}$

where ρ_mem is the dissociation polymer density; EW is the isomeric equivalents; ω is the volume fraction of the ionomers in the CL; D_nmw is the ionomer hydration water diffusion coefficient; λ_nf is the ionomer hydration water content.

(vii) Charge conservation

(37) $\nabla \left({\kappa }_{ele}^{eff}\nabla {\varphi }_{ele}\right)+S_{ele}=0$

(38) $\nabla \left({\kappa }_{ion}^{eff}\nabla {\varphi }_{ion}\right)+S_{ion}=0$

(39) ${\kappa }_{ele}^{eff}={\left(1-\epsilon -\omega \right)}^{1.5}{\kappa }_{ele}$

(40) ${\kappa }_{ion}^{eff}={\omega }^{1.5}{\kappa }_{ion}$

(41) ${\kappa }_{ion}=\left(0.5139{\lambda }_{nf}-0.326\right)\exp \left[2222\left(\frac{1}{303.15}-\frac{1}{T}\right)\right]$

where Κ_ion is the ionic conductivity; Κ_ele is the electron conductivity; ${\kappa }_{ion}^{eff}$ is the reference ion conductivity; ${\kappa }_{ele}^{eff}$ is the reference to the electron conductivity.

(viii) Energy conservation

(42) $\frac{\partial }{\partial t}\left({\left(\rho C_p\right)}_{fl,sl}^{eff}T\right)+\nabla {\left(\rho C_p\right)}_{fl}^{eff}{\overrightarrow{u}}_{g}T-\nabla {\kappa }_{fl.sl}^{eff}\nabla T=S_T$

(43) ${\left(\rho C_p\right)}_{fl,sl}^{eff}=\epsilon \left[s_l{\rho }_l{\left(C_p\right)}_l+\left(1-s_l\right){\rho }_g{\left(C_p\right)}_g\right]+\left(1-\epsilon -\omega \right){\rho }_{sld}{\left(C_p\right)}_{sld}+\omega {\rho }_{mem}{\left(C_p\right)}_{mem}$

(44) ${\left(\rho C_p\right)}_{fl}^{eff}$

(45) ${\left(C_p\right)}_g$

(46) $k_{fl,sl}^{eff}=\epsilon \left[s_l{k}_l+\left(1-s_l\right)k_g\right]+\left(1-\epsilon -\omega \right)k_{sld}+\omega k_{mem}$

(47) $k_g={\displaystyle \sum _i\frac{X_i{k}_i}{{\displaystyle {\sum }_j{X}_i{\psi }_{kij}}}}$

(48) ${\psi }_{kij}=\frac{{\left[1+{\left(\frac{k_i}{k_j}\right)}^{0.5}{\left(\frac{M_i}{M_j}\right)}^{0.25}\right]}^{0.5}}{{\left[8\left(1+\frac{M_i}{M_j}\right)\right]}^{0.5}}$

where ρ_sld is the solid phase density (solid phase excluding membrane electrolytes and ice); (C_p)_sld is the effective volume heat capacity (solid phase excluding membrane electrolytes); ${\left(\rho C_p\right)}_{fl,sl}^{eff}$ is the effective volume heat capacity (solid phase including membrane electrolytes); ${\left(\rho C_p\right)}_{fl}^{eff}$ is the effective volume heat capacity (liquid phase); (C_p)_g is the specific heat capacity of the mixed gas; (C_p)_l is the liquid water specific heat; (C_p)_ice is the ice specific heat; (C_p)_mem is the ionomer specific heat; (C_p)_sld is the solid phase (not include the membrane electrolytes and the ice) specific heat; $k_{fl,sl}^{eff}$ is the effective thermal conductivity; k_l is the thermal conductivity of the liquid water; k_g is the thermal conductivity of the mixed gas; k_sld is the solid phase does not include the membrane electrolyte thermal conductivity; k_mem is the thermal conductivity of the ionomers.

Equations of phase change

In the actual work of the fuel cell, the mutual phase transition between liquid water, water vapor, and ionomer hydration water is involved. Generally, it is expressed by the following equations:

(i) Water vapor-liquid water

(49) $S_{v-l}=\{\begin{array}{ccc}{\gamma }_{cond}\epsilon \left(1-s_l\right)\frac{\left(p_{vp}-p_{sat}\right)M_{H_{2}O}}{RT}& \text{if}& p_{vp}\ge p_{sat}\\ {\gamma }_{evap}\epsilon s_l\frac{\left(p_{vp}-p_{sat}\right)M_{H_{2}O}}{RT}& \text{if}& p_{vp}<p_{sat}\end{array}$

(50) ${\log }_{10}\left(\frac{p_{sat}}{101325}\right)=-2.1794+0.02953\left(T-273.15\right)-9.1837\times {10}^{-5}{\left(T-273.15\right)}^2+1.4454\times {10}^{-7}{\left(T-273.15\right)}^3$

where S_v–l is the steam-liquid water source item; λ_evap is the condensing rate; λ_cond is the condensing rate.

(ii) Ionomer hydration water-water vapor

(51) $S_{n-v}=\{\begin{array}{c}{\zeta }_{n-v}\frac{{\rho }_{mem}}{EW}\left({\lambda }_{nf}-{\lambda }_{eq}\right)\left(1-s_l\right)\text{if}{\lambda }_{nf}\ge {\lambda }_{eq}\\ {\zeta }_{n-v}\frac{{\rho }_{mem}}{EW}\left({\lambda }_{eq}-{\lambda }_{nf}\right)\left(1-s_l\right)\text{if}{\lambda }_{nf}<{\lambda }_{eq}\end{array}$

(52) ${\lambda }_{eq}=\{\begin{array}{ccc}0.043+17.81a-39.85a^2+36.0a^3& \text{if}& 0\le a\le 1\\ 14.0+1.4\left(a-1\right)& & 1<a\le 3\end{array}$

(53) $a=\frac{p_{vp}}{p_{sat}}+2S_l$

where S_n-v is the ionomer hydration water-water vapor source term; ζ_n–v is the ionomer hydration water-water vapor conversion rate; ζ_v–n is the water vapor-ionomer hydration water conversion rate; p_sat is the saturated vapor pressure.

Source phase

The individual source phases in the formula are calculated by the following equations:

(i) Material source item

(54) $S_m=S_{H_2}+S_{O_2}+S_{vp}$

(55) $S_{H_2}=\{\begin{array}{ccc}-\frac{j_a}{2F}{M}_{H_2}& \text{in}& \begin{array}{cc}\text{anode}& \text{CL}\end{array}\\ 0& & \begin{array}{cc}\text{other}& \end{array}\end{array}$

(56) $S_{O_2}=\{\begin{array}{ccc}-\frac{j_c}{4F}{M}_{O_2}& \text{in}& \begin{array}{cc}\text{cathode}& \text{CL}\end{array}\\ 0& & \begin{array}{cc}\text{other}& \end{array}\end{array}$

(57) $S_{vp}=\{\begin{array}{cc}-S_{v-l}+S_{n-v}{M}_{H_{2}O}& \text{CL}\\ -S_{v-l}& \text{other}\end{array}$

(ii) Momentum source term

(58) $S_u=\{\begin{array}{ccccc}-\frac{{\mu }_g}{K_g}{\overrightarrow{u}}_g& \text{in}& \text{CL}& \text{and}& \text{GDL}\\ 0& & & \text{other}& \end{array}$

(iii) Liquid water source item

(59) $S_l=S_{v-l}$

(iv) Ionomer hydration water source item

(60) $S_{nmw}=\{\begin{array}{cc}0& \text{membrane}\\ \frac{j_c}{2F}-S_{n-v}+S_{EOD}& \begin{array}{cc}\text{cathode}& \text{CL}\end{array}\\ -S_{n-v}+S_{EOD}& \begin{array}{cc}\text{anode}& \text{CL}\end{array}\end{array}$

(61) $S_{EOD}=\nabla \left(\frac{n_d}{F}{\kappa }_{ion}^{eff}\nabla {\varphi }_{ion}\right)$

(62) $n_d=\frac{2.5{\lambda }_{nf}}{22}$

where n_d is the electro-seepage drag coefficient; S_EOD is the drag source phase.

(v) Electronic source item

(63) $S_{ele}=\{\begin{array}{ccc}-j_a& \text{anode}& \text{CL}\\ j_c& \text{cathode}& \text{CL}\\ 0& \text{other}& \end{array}$

where j_a is the anodic volume-exchange current density; j_c is the cathodic volume-exchange current density.

(vi) Proton source term

(64) $S_{ion}=\{\begin{array}{ccc}{j}_a& \text{anode}& \text{CL}\\ -j_c& \text{cathode}& \text{CL}\\ 0& \text{other}& \end{array}$

(65) $j_a=\left(1-s_l\right)j_{0,a}^{ref}{\left(\frac{c_{H_2}}{c_{H_2}^{ref}}\right)}^{0.5}\times \left[\exp \left(\frac{2{\alpha }_{a}F}{RT}{\eta }_{act,a}\right)-\exp \left(-\frac{2{\alpha }_{c}F}{RT}{\eta }_{act,a}\right)\right]$

(66) $j_c=\left(1-s_l\right)j_{0,c}^{ref}\left(\frac{c_{O_2}}{c_{O_2}^{ref}}\right)\times \left[-\exp \left(\frac{4{\alpha }_{a}F}{RT}{\eta }_{act,c}\right)+\exp \left(-\frac{2{\alpha }_{c}F}{RT}{\eta }_{act,c}\right)\right]$

(67) $j_{0,a}^{ref}={10}^9\exp \left[-1400\left(\frac{1}{T}-\frac{1}{353.15}\right)\right]$

(68) $j_{0,c}^{ref}={10}^4\exp \left[-7900\left(\frac{1}{T}-\frac{1}{353.15}\right)\right]$

where $j_{0,a}^{ref}$ is the anode volume reference exchange current density; $j_{0,c}^{ref}$ is the cathodic volume reference exchange current density.

(vii) Energy source term

(69) $S_T=\{\begin{array}{cc}{j}_a\left|{\eta }_{act}\right|+{\Vert \nabla {\varphi }_{ele}\Vert }^2{\kappa }_{ele}^{eff}+{\Vert \nabla {\varphi }_{ion}\Vert }^2{\kappa }_{ion}^{eff}+S_{pc}& \text{anode}\text{CL}\\ -\frac{j_{c}T\Delta S}{2F}+j_c\left|{\eta }_{act}\right|+{\Vert \nabla {\varphi }_{ele}\Vert }^2{\kappa }_{ele}^{eff}+{\Vert \nabla {\varphi }_{ion}\Vert }^2{\kappa }_{ion}^{eff}+S_{pc}& \text{cathode}\text{CL}\\ {\Vert \nabla {\varphi }_{ele}\Vert }^2{\kappa }_{ele}^{eff}+S_{pc}& \text{GDL}\left(69\right)\\ {\Vert \nabla {\varphi }_{ele}\Vert }^2{\kappa }_{ele}^{eff}& \text{BP}\\ {\Vert \nabla {\varphi }_{ion}\Vert }^2{\kappa }_{ion}^{eff}+S_{pc}& \text{membrane}\\ 0& \text{other}\end{array}$

(70) $S_{pc}=\{\begin{array}{cc}{h}_{cond}{S}_{v-l}& \text{GDL}\\ h_{cond}\left(S_{v-l}-S_{n-v}{M}_{H_{2}O}\right)& \text{CL}\\ 0& \text{other}\end{array}$

where η_act is the overpotential; S_pc is the latent heat function; h_cond is the latent heat of condensation.

Computational domain

The model in this paper is a single cell with a 3 cm × 3 cm flow field. The calculation area includes the proton exchange membrane, cathode catalytic layer (CL), anode CL, cathode GDL, anode GDL, cathode flow field, and anode flow field. The original flow field dimension is width × height = 1 mm × 1 mm (Figure 1a). The precision dimensions are 0.7 mm × 0.7 mm (Figure 1b) and 0.5 mm × 0.5 mm (Figure 1c). In addition, the mechanism of water and heat transfer within the internal cross-section of the flow field is also discussed in detail, and the cross-section is shown in Figure 1d.

Figure 1
The calculated area

Model hypothesis

• (i) The fuel cell reaction process is a constant temperature process;

• (ii) Flow is laminar, an incompressible flow;

• (iii) The gravity effect is ignored;

• (vi) Porous media, such as GDL and CL, are all isotropic;

• (v) GDL and CL are hydrophobic, and their contact angle is greater than 90°;

• (vi) When the liquid water crosses the junction of the GDL, the gas carries the flow field, that is, liquid water only exists in GDL and CL.

Basic parameters

In this paper, the simulation model was built by the COMSOL Multiphysics 6.0, and the parameters of the model are shown in Table 1.

Table 1
Parameters

Parameter	Value
Cell voltage	1 V
Reference pressure	1.0133 × 105 Pa
Anode equilibrium potential	0.030549 V
Anode power viscosity	2.46 × 10-5 Pa s
Anode stoichiometry ratio	2.5
Anode transfer coefficient	0.5
Cathode equilibrium potential	1.2069 V
Cathode power viscosity	1.19 × 10-5 Pa s
Cathode stoichiometric ratio	1.5
Cathode transfer coefficient	0.5
Cell length	2 mm
Channel height	0.05 mm
Channel width	0.5 mm
CL conductivity	25 S m-1
CL porosity	0.3
CL thickness	18 µm
Density of liquid water	1000 kg m-3
Diffusivity of hydrogen in the ionomers	2 × 10-9 m2 s-1
Dissociation polymer density	2000 kg m-3
Electrode volume fractions, CL	0.4
Electrolyte liquid phase volume fraction	0.3
Exchange current density, hydrogen oxidation	100 A m-2
Exchange current density, oxygen reduction	0.01 A m-2
Exchange membrane with wet molar volume	1.8 × 10-4 m3 mol 1
Exchange the membrane for a dry molar volume	5.5 × 10-4 m3 mol 1
Gas constant	8.314 J mol-1 K-1
GDL conductivity	100 S m-1
GDL permeability	6.2 × 10-12 m2
GDL porosity	0.5
GDL thickness	215 µm
Hydrogen gas reference diffusion coefficient	1.24 × 10-4 m2 s1
Hydrogen molar mass	2 g mol-1
Hydrogen-water vapor binary diffusion coefficient	9.15 × 10-5 m2 s1
Initial humidity of the anode	100%
Initial humidity of the cathode	100%
Isomeric equivalents	1.1 kg mol-1
Liquid water density	1000 kg m-3
Liquid water viscosity	3.7 × 10-4 Pa s
Mixed gas viscosity	2.46 × 10-5 Pa s
Nitrogen molar mass	28 g mol-1
Nitrogen-water vapor binary diffusion coefficient	2.56 × 10-5 m2 s1
Number of integrals of the ionomer	0.2
Oxygen molar mass	32 g mol-1
Oxygen reference diffusion coefficient	2.8 × 10-5 m2 s-1
Oxygen-nitrogen gas binary diffusion coefficient	2.2 × 10-5 m2 s-1
Oxygen-water vapor binary-dependent diffusion coefficient	2.82 × 10-5 m2 s1
Porous electrode permeability	1.24 × 10-12 m2
Proton exchange membrane with an equivalent molar volume	5.17 × 10-4 m3 mol 1
Proton-exchange membrane thickness	15 µm
Reference diffusion coefficient of water in the anode	1.24 × 10-4 m2 s1
Reference diffusion coefficient of water in the cathode	3.6 × 105 m2 s-1
Reference temperature	353.15 K
Rib width	1 mm
Solid-phase volume fraction, GDL	0.4
Water molar mass	18 g mol-1
Water vapor viscosity	2.1 × 105 Pa s

• GDL: gas diffusion layer; CL: catalytic layer.

Boundary condition

Before the conserved equation could be solved, the boundary conditions and the initial values must be set. For the flow boundaries, the entrance mass flow rate is determined by the respective current density, gas stoichiometry, and active reaction area.

The total cathode-molar mass:

(71) $M_c=x_{O_{2,in}}\times M_{w,O_2}+x_{H_2{O}_{c,in}}\times M_{w,H_{2}O}+x_{N_{2,in}}\times M_{w,N_2}$

The total anode-molar mass:

(72) $M_a=x_{H_{2,in}}\times M_{w,H_2}+x_{H_2{O}_{a,in}}\times M_{w,H_{2}O}$

The cathode mass flow rate:

(73) $m_c=\frac{stoic{h}_c\times i\times A_{cl}\times M_c}{4F{x}_{O_{2,in}}}$

The anode molar flow:

(74) $m_a=\frac{stoic{h}_a\times i\times A_{cl}\times M_a}{2F{x}_{H_{2,in}}}$

where x_O2,in is the inlet mole fraction of oxygen; x_H2Oc,in is the inlet mole fraction of cathode water vapor; x_N2,in is the inlet mole fraction of nitrogen gas; M_w,O2 is the molar masses of oxygen; M_w,H2O is the molar mass of water; M_w,N2 is the molar mass of nitrogen; x_H2,in is the inlet molar fraction of hydrogen; x_H2Oc,in is the inlet molar fraction of anode water vapor; stoich_c is the cathode stoichiometric ratio; i is the current density; A_cl is the active reaction area; stoich_a is the anode stoichiometry ratio.

Furthermore, the inlet velocity was expressed in accordance with the channel inlet area and the calculated mass flow rate.

For the charge conservation equation, the potential of the anode GDL contacting the bipolar plate was set as the reference voltage (0 V). At the interface of the cathode GDL and the bipolar plate, the potential was defined as the cell operating voltage. Furthermore, the model used a constant temperature for all the calculated fields. For the remaining boundary conditions, no-flux or symmetry conditions were used as follows:

(75) $\frac{\partial \chi }{\partial n}=0$

where χ is the variable associated with the model solution.

Model verification

The experimental procedure in this study begins with the activation of the membrane electrode assembly (MEA), employing an intermittent activation approach with a duration of twelve hours each day for a total of four days. Following activation, the fuel cell performance testing is conducted by discharging the cell from 0.4 to 0.8 V in increments of 0.05 V. During testing, an assembly torque of 3.0 N m is applied to the proton exchange membrane fuel cell (PEMFC), and the experimental cell temperature is maintained at 70 °C.

Figure 2b shows the comparison diagram of the polarization curve simulated and experimental fitted under the same conditions and the same geometric parameters (the specific parameters are shown in Table 2) (the experimental setup diagram is shown in Figure 2a). Figure 2 shows that the simulated polarization curve is highly consistent with the experiment, hence the reliability of the simulation model used in this paper.

Figure 2
Test verification: (a) schematic diagram of experiment device; (b) comparison of experiment and simulation results

Table 2
Conditional parameters of the experimental PEMFC

Parameter / unit	Value
GDL porosity	0.4
GDL thickness / mm	0.215
CL thickness / mm	0.018
CL porosity	0.5
Platinum load / mg cm-2	0.5
Flow field type	Three snake-shaped flow field
Effective area / cm2	25
Air flow rate / (mL min-1)	680
Hydrogen flow rate / (mL min-1)	340
Cell temperature / °C	70
Gas temperature / °C	70

• GDL: gas diffusion layer; CL: catalytic layer.

RESULTS AND DISCUSSION

Electrochemical characteristics of different precision flow fields

The polarization curve is a key tool to evaluate and optimize the performance of fuel cells. Generally, the polarization curve can be divided into three regions: the activated polarization region of low current density, the ohmic polarization region of medium current density, and the concentration polarization region of high current density. The above three regions represent the electrochemical reaction characteristics, ohm impedance characteristics, and mass transfer characteristics of the fuel cell, respectively. As can be seen in Figure 3, the activation polarization of the three precision flow fields is almost the same because the flow field precision design presented in this paper does not change the catalytic efficiency of the catalyst required for the electrochemical reaction.

Figure 3
Comparison of PEMFC polarization curves of flow fields with different precision

With the improvement of flow precision (the smaller the flow size, the higher the precision), both ohm polarization and concentration polarization are reduced, and the improvement effect of concentration polarization is more obvious. The reduction of ohmic polarization shows that the precision design can effectively reduce the contact resistance between the flow channel and GDL and reduce the ohmic overpotential. The improvement of the concentration difference polarization effect shows that the precision design of the flow channel can improve the mass transfer performance and water removal characteristics, effectively improve the utilization rate of gas, and finally reduce the overpotential of the concentration difference.

Concentration overpotential is a polarization phenomenon of fuel cells, which usually refers to the voltage loss on the electrode surface in the cell due to the difference in material concentration. The concentration difference overpotential mainly involves the transmission process of electrolyte, fuel, and oxygen inside the cell. Figure 4 shows the concentration overpotential of different precision flow fields. Obviously, with the increase of cell voltage, the concentration overpotential increases significantly. In addition, with the improvement of the cell precision degree, the concentration difference overpotential is significantly improved. As can be seen in this the figure, when the flow channel is refned from 1 mm × 1mm to 0.7 mm × 0.7 mm, the concentration difference overpotential at 0.4 V decreases from 0.3554 to 0.1033, a significant decrease. When the precision degree is further increased to 0.5 mm × 0.5 mm, the concentration difference overpotential at 0.4 V was further increased, but the improvement degree decreased.

Figure 4
Concentration difference overpotential of different precision flow fields of PEMFC

The reason for the above phenomenon is that when the degree of precision is higher, the gas distribution of the flow channel is more uniform. As shown in Figure 1, with the improvement of precision, the number of channels increases, and the flow area through the GDL also increases. On the other hand, it can improve the gas distribution inside the GDL. The lower concentration difference overpotential can lift the higher material transmission and reaction rate of the fuel cell, which enables the cell to react under the condition of the lower voltage, thus improving the efficiency of the cell.

Heat and water transfer characteristics of different precision flow fields

Figure 5 shows the oxygen inhomogeneity C_v on the surface of PEMFC cathode GDL, which is usually used to represent the difference or imbalance between different parts. The calculation formula is:

Figure 5
Analysis of oxygen inhomogeneity on the surface of cathode GDL of different precision flow field PEMFC

(76) $\begin{array}{c}{C}_v=\frac{sd}{mean}\\ sd=\sqrt{{\displaystyle \sum \frac{{\left(x_i-\overline{x}\right)}^2}{n}}}\\ mean=\frac{{\displaystyle \sum x_i}}{n}\end{array}$

where sd is the standard deviation, mean is the average value.

As can be seen in Figure 5, with the increase of the cell voltage, the C_v of the cathode GDL surface continuously increases, and the slope also gradually increases. The improvement of C_v indicates that oxygen distribution’s inhomogeneity on the cathode GDL surface increases and greater inhomogeneity occurs in the low voltage area (high current density) area. Moreover, increasing flow channel precision improves the oxygen inhomogeneity on the cathode GDL surface.

The above phenomenon is because oxygen diffusion on the cathode surface is limited at a high current density. The cell needs more oxygen when the current density increases to support the oxygen reduction reaction. However, the oxygen transfer rate in the gas channel or electrolyte may not meet this demand, resulting in a lower oxygen concentration near the flow outlet, thus causing an uneven distribution. Therefore, the precision design of the flow channel can improve the diffusion characteristics of oxygen in the GDL, thus reducing the inhomogeneity of oxygen.

Figure 6 shows the molar fraction of oxygen on the surface of the PEMFC cathode GDL with different precision flow fields at the voltage of 0.4 V. Obviously, the oxygen distribution on the cathode GDL surface is more uniform with the increasing precision. Specifically, a significant region of low molar fraction of oxygen appeared at the GDL below the rib plate, especially in the 1 mm × 1 mm flow field PEMFC. This phenomenon is significantly improved with the improvement of flow field precision.

Figure 6
Molar fraction of oxygen on the cathode GDL of PEMFC with different precision flow fields (cell voltage = 0.4 V): (a) 1 mm × 1 mm; (b) 0.7 mm × 0.7 mm; (c) 0.5 mm × 0.5 mm

The reason for the above phenomenon is that the precision design of the flow channel can reduce the flow resistance of the gas in the channel and the pressure drop. In addition, the increased number of flow channels with high precision PEMFC can reduce the gas flow speed difference within the flow channel and prevent the gas flow rate in some regions from being much higher than others. Eventually, the uniformity of the gas distribution is improved.

Figure 7 shows the current density inhomogeneity of the PEMFC cathode CL surface with different precision. It can be seen that the current density inhomogeneity of the PEMFC cathode CL surface gradually decreases with the improvement of flow channel precision. The uneven distribution of current density will make the current density of some parts too high, leading to a low local electrochemical reaction rate, thus affecting the cell’s overall performance. Moreover, too high a local current density will cause local hot spots, leading to the acceleration of electrochemical reactions in the high region and the aging of this region. In addition, long-term local hotspots can lead to perforation of the proton exchange membrane, ultimately causing cell damage. Therefore, the inhomogeneity of the current density distribution of PEMFC can be improved, thus improving the cell life.

Figure 7
Current density inhomogeneity on the surface of the cathode CL of flow field

Figure 8 shows the surface current density distribution of PEMFC cathode CL at 0.4 V. It is obvious from the figure that the cathode CL surface current density distribution uniformity increases with increasing precision, and the maximum and minimum difference of current density of 1 mm × 1 mm, 0.7 mm × 0.7 mm, and 0.5 mm × 0.5 mm are 957, 724 and 558 mA cm², respectively. The precision design can improve the current density’s distribution uniformity. On the one hand, the uniformity of the gas distribution makes the gas reach the active area surface and makes the chemical reaction more uniform. On the other hand, the precision design of the flow field can improve the water removal performance of the flow channel, discharge the water produced by the reaction in time, and effectively prevent the local water flooding from hindering the mass transfer of the reaction gas.

Figure 8
Current density distribution of PEMFC cathode CL surfaces with different precision flow fields (cell voltage = 0.4 V): (a) 1 mm × 1 mm; (b) 0.7 mm × 0.7 mm; (c) 0.5 mm × 0.5 mm

Figure 9 shows the air mass transfer velocity of PEMFC at 0.4 V (vertical direction). The figure shows positive and negative mass transfer velocity because both the intake and exhaust (the gas carries liquid water out) processes exist in the GDL. The mass transfer speed of air in GDL is significantly accelerated with the increasing precision of the flow field. The reason is that the precision design of the channel can reduce the spacing between the channels and reduce the path of gas diffusion, thus improving the mass transmission speed. In addition, when the inlet flow rate is certain, decreasing the inlet channel’s surface area increases the gas’s flow rate, thus increasing its mass transfer speed in the GDL.

Figure 9
Mass transfer velocity of PEMFC with different precision flow fields (cell voltage = 0.4 V): (a) 1 mm × 1 mm; (b) 0.7 mm × 0.7 mm; (c) 0.5 mm × 0.5 mm

Analysis of PEMFC phase transition characteristics of different precision flow fields

Within the normal operating temperature range of the fuel cell (60-80 °C), there are three states of water: water vapor, liquid water, and ionomer hydration water. Under certain conditions, water transforms in the three states. For water vapor and liquid water, the phase transition between them is determined by the water vapor pressure and the saturated vapor pressure. According to Equations 49 and 50, when the water vapor pressure is greater than the saturated vapor pressure, the water vapor will change toward liquid water. On the contrary, if the water vapor pressure is less than the saturated vapor pressure, the liquid water will evaporate the phase change and then turn into water vapor.

For ionomer hydration water and liquid water, the phase transition between them is influenced by the relationship between the equilibrium water content and the ionomer hydration water content. According to Equation 51, water vapor and ionomer hydration water are interconverted due to the relationship between equilibrium water content and ionomer hydration water content. This section will focus on the influence of the flow field precision design on the water phase transition in PEMFC.

In fuel cells, liquid water is the main obstacle to gas mass transfer, especially in porous media; liquid water blocking the pore structure will cause a serious impact on the intake of air. Figure 10 shows the liquid water saturation of GDL in the adjacent flow field PEMFC. As can be seen in the figure, the saturation of liquid water inside GDL gradually decreases with the improvement of the flow precision. The reason is that, on the one hand, the flow field precision design can improve the mass transfer characteristics of the gas inside the porous medium, thus allowing the gas to take away more liquid water. On the other hand, the precision design of the flow channel can make the electrochemical reaction more sufficient, increase the consumption of oxygen, and reduce the pressure of oxygen, thus increasing the pressure of cathode water vapor and promoting the phase transition of liquid water to the direction of water vapor.

Figure 10
Distribution of GDL liquid water in different precision flow field PEMFC (equal scale amplification, same size of MEA in the three figures): (a) 1 mm × 1 mm; (b) 0.7 mm × 0.7 mm; (c) 0.5 mm × 0.5 mm

Figure 11 shows the hydration water distribution of ionomer domains between the adjacent flow channels in PEMFC fields. Obviously, with the increase in flow field precision, hydration water gradually decreased. The reason is that with increasing flow field precision, the saturation of the liquid water inside the porous medium decreases. According to Equations 52 and 53, with the decrease of liquid water saturation, the water activity a decreases, thus reducing the equilibrium water content and promoting the phase transition from ionomer hydration water to water vapor. With the vapor phase transfer of ionomer hydration water, vapor pressure gradually increases. Finally, the removal of the liquid water is promoted under evaporation.

Figure 11
Hydration water distribution of adjacent flow fields of the different precision flow field (equal scale amplification, same size of MEA in the three figures): (a) 1 mm × 1 mm; (b) 0.7 mm × 0.7 mm; (c) 0.5 mm × 0.5 mm

CONCLUSIONS

This paper establishes a three-dimensional multiphase model to explore the influence of flow field precision design on the electrochemical characteristics, water and heat transfer characteristics, and fuel cell phase change characteristics.

With the improvement of the flow field precision, the concentration difference overpotential is significantly improved, the gas utilization rate is improved, and the water removal performance of the flow field is increased. In this paper, the inhomogeneity (C_v) was used to measure the gas and current density distribution characteristics of different precision flow fields. With increasing precision, the oxygen distribution is more uniform, and the region of low molar fraction of oxygen is significantly reduced. Meanwhile, the current density distribution is more uniform, and the difference between maximum and minimum current density decreases. With the improvement of the flow field precision, the mass transfer speed of the gas is significantly improved, which increases the mass transfer characteristics of the gas and is conducive to removing liquid water. With the improvement of the flow field precision, the liquid water content inside GDL decreases, and the ionomer hydration water decreases, indicating that the flow field precision design can improve the water removal performance of the flow field.

REFERÊNCIAS

• ¹
  Yuan, H.; Dai, H.; Ming, P.; Zhao, L.; Tang, W.; Wei, X.; Chem. Eng. J. 2022, 431, 134035. [Crossref]
  » Crossref
• ²
  Tang, X.; Zhang, Y.; Xu, S.; Int. J. Heat Mass Transfer 2023, 206, 123966. [Crossref]
  » Crossref
• ³
  Shi, L.; Liu, P.; Zheng, M.; Xu, S.; Energy Convers. Manage. 2023, 278, 116725. [Crossref]
  » Crossref
• ⁴
  Cheng, C.; Yang, Z.; Liu, Z.; Tongsh, C.; Zhang, G.; Xie, B.; He, S.; Jiao, K.; Appl. Energy 2021, 302, 117555. [Crossref]
  » Crossref
• ⁵
  Khazaee, I.; Sabadbafan, H.; Heat Mass Transfer 2016, 52, 993. [Crossref]
  » Crossref
• ⁶
  Chowdhury, M. Z.; Genc, O.; Toros, S.; Int. J. Hydrogen Energy 2018, 43, 10798. [Crossref]
  » Crossref
• ⁷
  Ramin, F.; Sadeghifar, H.; Torkavannejad, A.; Int. J. Heat Mass Transfer 2019, 129, 1151. [Crossref]
  » Crossref
• ⁸
  Timurkutluk, B.; Chowdhury, M. Z.; Fuel Cells 2018, 18, 441. [Crossref]
  » Crossref
• ⁹
  Ashrafi, M.; Shams, M.; Appl. Energy 2017, 208, 1083. [Crossref]
  » Crossref
• ¹⁰
  Han, I. S.; Lim, J.; Jeong, J.; Shin, H. K.; Renewable Energy 2013, 54, 180. [Crossref]
  » Crossref
• ¹¹
  Niu, Z.; Fan, L.; Bao, Z.; Jiao, K.; Int. J. Energy Res. 2018, 42, 3328. [Crossref]
  » Crossref
• ¹²
  Zhang, G.; Xie, B.; Bao, Z.; Niu, Z.; Jiao, K.; Int. J. Energy Res. 2018, 42, 4697. [Crossref]
  » Crossref
• ¹³
  Ghadhban, S. A.; Alawee, W. H.; Dhahad, H. A.; Case Studies in Thermal Engineering 2021, 24, 100841. [Crossref]
  » Crossref
• ¹⁴
  Asadzade, M.; Shamloo, A.; Int. J. Energy Res. 2017, 41, 1730. [Crossref]
  » Crossref
• ¹⁵
  Kang, H. C.; Jum, K. M.; Sohn, Y. J.; Int. J. Hydrogen Energy 2019, 44, 24036. [Crossref]
  » Crossref
• ¹⁶
  Damian-Ascencio, C. E.; Saldana-Robles, A.; Hernandez-Guerrero, A.; Cano-Andrade, S.; Energy 2017, 133, 306. [Crossref]
  » Crossref
• ¹⁷
  Cano-Andrade, S.; Hernandez-Guerrero, A.; von Spakovsky, M. R.; Damian-Ascencio, C. E.; Rubio-Arana, J. C.; Energy 2010, 35, 920. [Crossref]
  » Crossref
• ¹⁸
  Geneve, T.; Regnier, J.; Turpin, C.; Int. J. Hydrogen Energy 2016, 41, 516. [Crossref]
  » Crossref
• ¹⁹
  Wang, L. P.; Zhang, L. H.; Jiang, J. P.; Appl. Mech. Mater. 2011, 44-47, 2404. [Crossref]
  » Crossref
• ²⁰
  Su, A.; Weng, F. B.; Chi, P. H.; Lu, S. M.; Jung, G. B.; Tu, C. H.; Ferng, Y. M.; Proc. Inst. Mech. Eng., Part A 2007, 221, 617. [Crossref]
  » Crossref
• ²¹
  Tan, Q.; Lei, H.; Liu, Z.; Int. J. Hydrogen Energy 2022, 47, 11975. [Crossref]
  » Crossref
• ²²
  Shimpalee, S.; Van Zee, J. W.; Int. J. Hydrogen Energy 2007, 32, 842. [Crossref]
  » Crossref
• ²³
  Chen, S.; Xia, Z.; Zhang, X.; Wu, Y.; Procedia Eng. 2019, 158, 1678. [Crossref]
  » Crossref
• ²⁴
  Freire, L. S.; Antolini, E.; Linardi, M.; Santiago, E. I.; Passos, R. R.; Int. J. Hydrogen Energy 2014, 39, 12052. [Crossref]
  » Crossref
• ²⁵
  Kumar, P. M.; Kolar, A. K.; Int. J. Hydrogen Energy 2010, 35, 671. [Crossref]
  » Crossref
• ²⁶
  Kumar, A.; Reddy, R. G.; J. Power Sources 2003, 113, 11. [Crossref]
  » Crossref
• ²⁷
  Fontana, E.; Mancusi, E.; da Silva, A.; Mariani, V. C.; de Souza, A. A. U.; Souza, S. M. A.; Int. J. Heat Mass Transfer 2011, 54, 4462. [Crossref]
  » Crossref
• ²⁸
  Korkischko, I.; Carmo, B. S.; Fonseca, F. C.; Fuel Cells 2017, 17, 809. [Crossref]
  » Crossref
• ²⁹
  Peng, Y.; Mahyari, H. M.; Moshfegh, A.; Javadzadegan, A.; Toghraie, D.; Shams, M.; Rostami, S.; Int. Commun. Heat Mass Transfer 2020, 115, 104638. [Crossref]
  » Crossref

Publication Dates

• Publication in this collection
  13 May 2024
• Date of issue
  2024

History

• Received
  24 Oct 2023
• Accepted
  23 Jan 2024
• Published
  18 Apr 2024