进口温湿度对质子交换膜燃料电池输出性能的影响

摘要

为研究阴阳极气体进口温度和相对湿度（relative humidity，RH）对质子交换膜燃料电池输出性能的影响，构建了一维非等温的两相模型，并采用团聚体子模型考虑催化层微观结构的影响。结果表明：阴阳极相对湿度相同时，干工况（RHa 50% / RHc 50%）下工作温度对燃料电池输出性能的影响有限；不同工作温度时获得的最高电流密度所对应的RH工况不同，90℃时可在高阳极湿度/低阴极湿度工况（RHa 90% / RHc 50%）下获得最高的电流密度，而70℃适合于在高湿度环境（RHa 90% / RHc90 %）下获得更好的输出性能。然后，基于粒子群算法，优化得出阳极湿度恒为90%、且在低阴极RH 0~50%时，温度越高，获得最大功率密度的阴极湿度越小。工作温度和阴极RH分别控制在90℃和14.6%时可获得最高的功率密度，约为0.88 W/cm²。

关键词

质子交换膜燃料电池; 进口温度; 相对湿度; 粒子群算法

质子交换膜燃料电池（proton exchange membrane fuel cell， PEMFC，以下简称“燃料电池”），作为一种将氢气中的化学能直接转化为电能的能量转换装置，因其具有较高的功率密度、较低的噪声、零排放等显著优势，被认为是未来最有前途的发电及动力装置之一，现已成为世界各国研究的热点^［1］。影响燃料电池输出性能的因素有很多，包括多孔介质层的结构及润湿性、催化层的聚合物和催化剂含量、膜的材料特性以及燃料电池工作条件等。在这些影响因素中，燃料电池的工作条件（如进气压力、进气相对湿度、进气温度、进气量等）对燃料电池内部气体、液态水、温度分布等有着直接的关系^［2-4］，进而影响燃料电池的输出性能。

为了进一步提高燃料电池功率密度，有不少学者针对燃料电池的不同操作条件对输出性能的影响作了研究。Akitomo等人^［5］采用极化测试和X射线同步技术研究了气体加压对高温条件下工作的燃料电池输出性能的影响。Pei等人^［6］发现通过调节入口相对湿度（relative humidity，RH）或电池温度可以降低电压的额外波动，使电压输出基本保持不变。Janicka等人^［7］为了满足膜的润湿需求以及避免局部水淹，采用动态电化学阻抗谱获得了满足最佳工作条件下的气体RH和电流负载。除了实验研究手段，数值模拟凭借较低的研究成本以及较高的效率，已成为研究燃料电池的有力工具。焦魁等人^［8］通过构建三维非等温多相模型，研究了不同工作电压、温度、和RH下的电流密度空间分布情况。结果发现最高电流密度区域随着RH的变化而改变，而工作温度对电流密度分布影响较小。许思传等人^［9］通过构建一维燃料电池气液两相流模型，研究了进气RH对燃料电池输出性能的影响。结果发现在小电流密度条件下，进气RH对电池输出电压影响较大，而在大电流密度条件下影响较小。为了获得不同工作参数对燃料电池输出性能和气体浓度分布的影响，陈会翠等人^［10］同样建立了三维燃料电池模型，并根据正交原理，对比分析了包括阴阳极进口RH、压力、气体计量比以及温度在内的7个操作参数在不同电流负载条件下对输出性能的敏感程度。

对于通过增大气体计量比以及进口压力的方法，尽管有利于液态水的去除以及反应气体分压的提高，但是不可避免会增加额外的泵功损失，不利于燃料电池净功率密度的提升。而通过调控进口气体温度和RH，同样可以直接影响燃料电池的性能输出。由于饱和蒸气压与温度呈指数关系，高温低湿环境下燃料电池更易遭遇“膜干”现象，影响膜的质子导电性。相反，低温高湿又极易产生“水淹”，影响反应气体的扩散。目前，围绕同时调节进口气体的温湿度来获得较高的输出性能的研究还比较少，且大部分研究仍把工作温度限定在60~80 ℃。鉴于高温对提高电化学反应速率、缩小冷却系统尺寸、增强对CO的适应性等方面具有明显优势^［11］，因此有必要研究更高温度下，进口温湿度变化对燃料电池输出性能的影响。本文将通过构建一维、非等温、两相模型，并耦合催化层团聚体子模型来考虑催化层微观结构的影响。重点对比研究不同进口温湿度条件下，燃料电池的输出性能变化。同时采用粒子群算法，寻找能够获得最佳功率密度的温湿度条件，从而为燃料电池输出性能的提高，提供较合理的温湿度调节方案。

1 计算模型的构建与验证

1.1　计算对象

计算区域如图1所示，包括质子交换膜（GORE-SELECT® 735.18）、阴阳极催化层（catalyst layer，CL）、微孔层（micro porous layer，MPL）以及气体扩散层（gas diffusion layer，GDL，Toray H-060）。由于实际燃料电池每个部件的厚度都小于其长度和宽度，因此本文将模型简化为垂直膜方向（x轴）上的一维问题。各个部件厚度所对应的x轴坐标如图1所示。

图1  计算对象及团聚体结构

Fig.1  Computational domain and agglomerate structure

1.2　模型建立的方法

模型涉及的基本控制方程包括：质量守恒方程、动量守恒方程、组分输运方程、能量守恒方程、液态水输运方程、膜态水输运方程以及电荷守恒方程，如表1所示。在表1中： ρ、u、C、µ分别为密度、速度、质量浓度和粘度；M为物质的分子质量；i和J为电流密度和体积电流密度；s和λ分别代表液态水饱和度和膜态水量；S为守恒方程源项；n_d和D则为电渗拖曳系数以及各物质的扩散系数。计算域内的气液两相流动采用双流体方法描述，即分开求解液相和气相组分。动量守恒方程简化为广义达西方程形式。膜态水输运方程充分考虑了膜态水在阴阳极两侧的跨膜输运，包括电渗拖曳（electro-osmotic drag）、反扩散（back diffusion）、热渗透（thermo-osmosis drag）^［12］。值得指出的是，由于催化层和微孔层的纳米级孔隙结构给预测气体扩散过程带来了更大的难度，因此本模型结合克努森扩散（Knudsen diffusion）对催化层和微孔层内的气体扩散过程进行了修正^［13-14］。模型涉及的其他本构关系和传递关联式以及众多符号参数具体详见已发表的研究内容^［15］。为了更好地求解上述方程，做出了如下假设：气体为理想气体，在燃料电池内部为层流运动；电子和气体组分无法通过膜进行渗透；通过电化学反应产生的水为液态水；所有问题均在稳态下进行求解。上述控制方程基于有限差分法，采用MATLAB （R2018a）BVP函数进行求解，计算误差选为10^-8，以确保计算精度。

表1  控制方程及源项

Tab.1  Governing functions and source terms

名称	表达式	源项
质量守恒方程	$\nabla ·\left({\rho }_g{\stackrel{⃑}{u}}_g\right)=S_m$	$S_m=S_{{\mathrm{H}}_2}+S_{{\mathrm{O}}_2}+S_{{\mathrm{H}}_2\mathrm{O}}$
组分输运方程	$\nabla ·\left(C_g{x}_i{\stackrel{⃑}{u}}_g\right)=\nabla ·\left(C_g{D}_i^{\mathrm{e}\mathrm{f}\mathrm{f}}\nabla x_i\right)+S_i$	$S_{{\mathrm{H}}_2}=\left\{\begin{array}{c}-\frac{i}{2F}\mathrm{a}\mathrm{C}\mathrm{L}\\ 0\mathrm{o}\mathrm{t}\mathrm{h}\mathrm{e}\mathrm{r}\mathrm{s}\end{array}\right.S_{O_2}=\left\{\begin{array}{c}-\frac{i}{4F}\mathrm{c}\mathrm{C}\mathrm{L}\\ 0\mathrm{o}\mathrm{t}\mathrm{h}\mathrm{e}\mathrm{r}\mathrm{s}\end{array}\right.$ $S_{{\mathrm{H}}_2\mathrm{O}}=\left\{\begin{array}{c}-S_{vl}\mathrm{G}\mathrm{D}\mathrm{L}\mathrm{s}+\mathrm{M}\mathrm{P}\mathrm{L}\mathrm{s}\\ -S_{vl}-S_{vd}\mathrm{C}\mathrm{L}\mathrm{s}\end{array}\right.$
能量守恒方程	$\nabla ·\left(\epsilon s{\rho }_l{c}_{p,l}{\stackrel{⃑}{u}}_{l}T+\epsilon \left(1-s\right){\rho }_g{c}_{p,g}{\stackrel{⃑}{u}}_{g}T\right)=\nabla ·\left(-k\nabla T\right)+S_T$	$S_T=\left\{\begin{array}{cc}{H}_{ec}{S}_{vl}+{\sigma }_e{‖\nabla {\varphi }_e‖}^2& \mathrm{G}\mathrm{D}\mathrm{L}\mathrm{s}+\mathrm{M}\mathrm{P}\mathrm{L}\mathrm{s}\\ H_{ec}\left(S_{vl}+S_{vd}\right)+{\sigma }_e{‖\nabla {\varphi }_e‖}^2+{\sigma }_p{‖\nabla {\varphi }_p‖}^2+J_a\left|{\eta }_a\right|-\frac{J_a∆S_{\mathrm{H}\mathrm{O}\mathrm{R}}T}{2F}& \mathrm{a}\mathrm{C}\mathrm{L}\\ H_{ec}\left(S_{vl}+S_{vd}\right)+{\sigma }_e{‖\nabla {\varphi }_e‖}^2+{\sigma }_p{‖\nabla {\varphi }_p‖}^2+J_c\left|{\eta }_c\right|-\frac{J_c∆S_{\mathrm{O}\mathrm{R}\mathrm{R}}T}{2F}& \mathrm{c}\mathrm{C}\mathrm{L}\\ {\sigma }_p{‖\nabla {\varphi }_p‖}^2& \mathrm{M}\mathrm{E}\mathrm{M}\end{array}\right.$
液态水守恒方程	$-\nabla ·\left(\frac{{\rho }_l}{M_{{\mathrm{H}}_2\mathrm{O}}}\frac{K_l}{{\mu }_l}\frac{\mathrm{d}{P}_c}{\mathrm{d}s}\nabla s\right)=S_s$	$S_{vl}=\left\{\begin{array}{c}{\gamma }_{\mathrm{c}\mathrm{o}\mathrm{n}}\epsilon \left(1-s\right)\left(C_V-C_{\mathrm{s}\mathrm{a}\mathrm{t}}\right)C_V\ge C_{\mathrm{s}\mathrm{a}\mathrm{t}}\\ {\gamma }_{\mathrm{e}\mathrm{v}\mathrm{a}}\epsilon s\left(C_V-C_{\mathrm{s}\mathrm{a}\mathrm{t}}\right)C_V<C_{\mathrm{s}\mathrm{a}\mathrm{t}}\end{array}\right.$
膜态水守恒方程	$\nabla ·\left(\frac{i}{F}{n}_d\right)+\nabla ·\left(D_T\nabla T\right)=\nabla ·\left(\frac{{\rho }_m}{EW}{D}_{\lambda }\nabla \lambda \right)+S_{\lambda }$	$S_{\lambda }={\gamma }_d\frac{{\rho }_m}{M_m}\left({\lambda }_{eq}-\lambda \right)$
电荷守恒方程	$\left\{\begin{array}{c}\nabla ·\left({\sigma }_e\nabla {\varphi }_e\right)=S_e\mathrm{G}\mathrm{D}\mathrm{L}\mathrm{s}+\mathrm{M}\mathrm{P}\mathrm{L}\mathrm{s}+\mathrm{C}\mathrm{L}\mathrm{s}\\ \nabla ·\left({\sigma }_p\nabla {\varphi }_p\right)=S_p\mathrm{C}\mathrm{L}\mathrm{s}+\mathrm{M}\mathrm{E}\mathrm{M}\end{array}\right.$	$S_p=\left\{\begin{array}{c}{J}_a\mathrm{a}\mathrm{C}\mathrm{L}\\ -J_c\mathrm{c}\mathrm{C}\mathrm{L}\end{array}\right.$ $S_e=\left\{\begin{array}{c}-J_a\mathrm{a}\mathrm{C}\mathrm{L}\\ J_c\mathrm{c}\mathrm{C}\mathrm{L}\end{array}\right.$
达西方程	${\stackrel{⃑}{u}}_g=-\frac{K_g}{{\mu }_g}\nabla P_g$ ${\stackrel{⃑}{u}}_l=-\frac{K_l}{{\mu }_l}\nabla P_l$

此外，为了更好考虑宏观尺度下催化层的真实结构，本文采用球形团聚体子模型来计算体积电流密度。如图1所示，假设球形团聚体均匀分布在催化层内，并分别由聚合物膜以及液膜包裹。由于氧气扩散速率及还原反应速率相比于氢气较慢，因此球形团聚体子模型只在阴极催化层内考虑，具体方程详见文献［16］。对于边界条件的设置，分别对阴阳极流道与GDL的交接面施以相应的边界条件（工作温度$T_c=T_a$；电压${\varphi }_{e,c}=U$，${\varphi }_{e,a}=\mathrm{0}$）。根据前人的研究，液态水饱和度为$s_a=s_a=\mathrm{0.12}$。阴阳极两侧流道与GDL交接面对应的氢气与氧气的摩尔分数分别由式（1）、式（2）给出。

$$x_{{\mathrm{H}}_{\mathrm{2}},a}=\mathrm{1}-\frac{P_{\mathrm{s}\mathrm{a}\mathrm{t}}R{H}_a}{P_a}$$

（1）

$$x_{{\mathrm{O}}_{\mathrm{2}},c}=\mathrm{0.21}\left(\mathrm{1}-\frac{P_{\mathrm{s}\mathrm{a}\mathrm{t}}R{H}_c}{P_c}\right)$$

（2）

式中：$x_{{\mathrm{H}}_{\mathrm{2}},a}$和$x_{{\mathrm{O}}_{\mathrm{2}},c}$分别为氢气和氧气的摩尔分数；$R{H}_a$、$R{H}_c$、$P_a$和$P_c$分别代表阳极和阴极侧的进口相对湿度和混合气体压力；$P_{\mathrm{s}\mathrm{a}\mathrm{t}}$则为水蒸气饱和压力。

1.3　模型的验证

为了验证模型的有效性，本文从燃料电池输出性能与欧姆阻抗两方面对模型进行验证。图2给出了工作温度为90℃和70℃时的极化曲线对比结果，而90℃时的欧姆阻抗对比结果由表2给出。实验所得极化曲线和欧姆阻抗分别由G20 Greenlight Innovation燃料电池测试平台以及Gamry Interface 5000E电化学工作站测取。实验工况为阴阳极相对湿度各取50%（RHa/RHc： 50%），工作温度为90℃和70℃ （Ta/Tc： 90℃及70℃），进口压力恒为1.5 atm（Pa/Pc： 1.5 atm；1 atm=101 325 Pa）。从中可以看出，模拟所得燃料电池输出性能与实验测试结果对比较好。在该两个工况下，由在线实验得到的欧姆阻抗与通过模拟计算所得到的结果，二者的最大相对误差约9%。因此，本文所构建的模型是可靠的。

（a）  RHa/RHc: 50%, Ta/Tc: 90℃

（b）  RHa/RHc: 50%, Ta/Tc: 70℃

图2  不同工作温度实验与数值模拟结果对比

Fig.2  Comparison of experimental and numerical results at different temperatures

表2  90℃时的欧姆阻抗结果对比

Tab.2  Comparison in ohmic resistances at 90℃

电流密度/ （A/cm2）	实验结果 阻抗/（Ω·cm2）	模拟结果 阻抗/（Ω·cm2）
0.2	0.0661	0.0629
0.4	0.0682	0.0634
0.6	0.0713	0.0643
0.8	0.0697	0.0651

2 结果与讨论

2.1　工作温度对输出性能的影响

本文首先围绕阴阳极两侧相同RH条件下，不同工作温度对燃料电池输出性能的影响进行研究，如图3所示。由于燃料电池最高功率密度通常在工作电压为0.5~0.6 V之间，本文所涉及的工况也符合该现象，因此主要选取工作电压为0.5 V和0.6 V时对应的输出性能进行研究。由图3可以看出，湿工况下（RHa90%/RHc90%），工作温度越高，对应的电流密度越小。需要强调的是，工作温度低于60℃时，燃料电池输出性能同样降低。而在干工况下（RHa50%/RHc50%），不同工作温度对应的电流密度相差不大。当工作电压为0.6 V时，对应的电流密度在1.0 A/cm²上下浮动，工作温度约为70-80℃时的电流密度略微高于其他温度。

图3  相同RH和电压下工作温度与电流密度的关系

Fig. 3  Operating temperature as function of current density at the same RH and voltage

相同的RH下（见图4），工作温度与阴极CL/MPL交接面上的氧气摩尔分数可做出如下解释：当燃料电池的工作电压为0.5或0.6 V时，其输出性能主要受质子导电能力、反应位点附近氧气含量以及电化学反应活性的共同影响；由于温度与饱和水蒸气压力呈指数关系，温度越高，相同RH下对反应气体的稀释作用就越显著，出现了图4中工作温度越高，氧气摩尔分数越小的现象；尤其湿工况下，当工作温度从60℃提高至90℃时，CL/MPL交接面上的氧气摩尔分数大幅下降，这也是导致图3中电流密度在湿工况下由于温度升高造成其降低幅度（气体饥饿）高于干工况的主要原因。值得一提的是，在干工况下燃料电池输出性能主要受限于质子传导的能力，因此在所研究的工作电压下，不同工作温度对应的电流密度相差不大。考虑到温度升高有利于提高反应活性，对于稀释作用较低的干工况，最大电流密度主要位于工作温度在70~80℃的范围内。

图4  相同RH和电压下工作温度与阴极CL/MPL交接面氧气摩尔分数的关系

Fig.4  O₂ mole fraction of CL/MPL interface at various operating temperatures with the same RH and voltage

2.2　相对湿度对输出性能的影响

气体进口RH与燃料电池的质子导电性、水的状态和分布以及反应气体含量有着直接关系，因此研究相同温度下，不同阴阳极RH对燃料电池输出性能的影响很有必要。图5分别给出了工作温度为70℃和90℃时，四种RH工况下对应的燃料电池电流密度。不同工作温度下，不同的阴阳极RH工况对燃料电池的输出性能展现出了不一样的影响程度。对于高温工况，由于高阳极RH及低阴极RH（RHa90%/RHc50%）工况既能满足阳极到阴极的电渗拖曳需求，又能最大程度降低高温带来的对氧气的稀释作用，因此当电压为0.5和0.6 V时，最高电流密度均出现在该工况，分别可达1.54 A/cm²和1.32 A/cm²。而对于70℃工况，则在湿工况（RHa90%/RHc90%）下获得了最高电流密度，分别为1.72 A/cm²（0.5 V）和1.49 A/cm²（0.6 V）。

图5  相同电压和工作温度下不同阴阳极RH与电流密度的关系

Fig.5  Current density at various RH cases with the same operating temperature and voltage

由于90℃对应的水蒸气饱和压力（70.13 kPa）是70℃（31.2 kPa）时的将近2.5倍，因此二者对于湿度变化的敏感度有所差异。结合图6进一步分析可知，以湿工况（RHa90%/RHc90%）和干工况（RHa50%/RHc50%）在CL/MPL交接面的氧气摩尔分数为例，90℃时干工况对应的氧气摩尔分数约为0.145，明显高于湿工况的对应值（0.105）。而当工作温度为70℃时，干工况（0.175）与湿工况（0.15）之间的变化幅度仅为0.025。因此，工作温度越高，对湿度变化引起的反应气体浓度变化越敏感。相比于90℃工况，燃料电池更倾向于在70℃且高湿度的环境下获得较好的输出性能。

图6  不同阴阳极RH与阴极CL/MPL交接面氧气摩尔分数的关系

Fig.6  O₂ mole fraction of CL/MPL interface in various RH cases

2.3　低阴极湿度下的性能优化

为了维持质子交换膜的润湿性，以满足高效的质子传导需求，目前主要采用对进口气体进行外部加湿的方式。但是外部加湿系统增加了整个燃料电池系统的成本，降低了整体的能源效率，不利于燃料电池的商业化。随着车用燃料电池对高功率密度的更大追求，对供气进行外部加湿并不是最可行的方案。因此，针对低湿度条件下通过优化进气温湿度以获得燃料电池的最佳输出性能是有意义的。根据粒子群算法精度高、收敛快的特点^［17］，本文基于该算法对不同工作温度（70~90℃）及低阴极RH（RHc：0~50%）条件下的最佳功率密度进行优化。计算过程中，设置迭代次数为500次，粒子维数为2，每次迭代20个粒子，阳极RH恒为90%，以满足电渗拖曳的需求，具体优化目标和条件如下：

优化目标： $P$（$P=Ui$）

优化对象： $\mathrm{0.5}\mathrm{V}<U<\mathrm{0.6}\mathrm{V}$

$$\mathrm{0}<\mathrm{R}\mathrm{H}\mathrm{c}<\mathrm{50}\mathrm{\%}$$

表3给出了不同工作温度下获得最佳功率密度的阴极RH。在所研究的温度范围内，最佳功率密度对应的工作电压基本稳定在0.54~0.55 V之间，并且随着温度降低，该电压值有逐渐减小的趋势。相反，对应的阴极RH则逐渐增大。可以明显看出，当阳极RH为90%，且工作温度和阴极RH分别控制在90℃和14.6%时，可获得最高的功率密度，约为0.88 W/cm²。

表3  不同工作温度下获得最佳功率密度的阴极RH优化结果

Tab.3  Cathode RH optimization results for optimum power density at different operating temperatures

温度/℃	工作电压/V	阴极RH/%	最佳功率密度/（W/cm2）
90	0.546	14.6	0.877
85	0.545	17.3	0.865
80	0.543	19.4	0.853
75	0.541	22.1	0.840
70	0.539	26.4	0.827

结合图7给出的不同工作温度下的欧姆阻抗和CL/MPL交接面的氧气摩尔分数，进一步分析可知，尽管70℃工况下对应的氧气摩尔分数略高于90℃，但由于同样阳极RH下，后者在阳极侧的水蒸气分压显然更大，其对应的欧姆阻抗约为0.043 Ω·cm²，比70℃工况对应的欧姆阻抗还要小近0.1Ω·cm²。得益于阳极侧较高的水蒸气含量，使得电渗拖曳过程在高温工况下能够得到较好的维持，因此，在氧气含量相差不大的情况下，高温显得更有优势。

图7  不同工作温度下最佳功率密度对应的欧姆阻抗和阴极CL/MPL交接面氧气摩尔分数

Fig.7  Ohmic resistance and O₂ mole fraction of CL/MPL interface corresponding to the optimal power density at different operating temperatures

3 结论

本文针对进口气体温湿度对质子交换膜燃料电池输出性能的影响特性，构建了一维、非等温、两相数值模型进行了模拟研究。得出了以下结论：

（1） 阴阳极相对湿度相同时，高温不利于燃料电池在湿工况（RHa 90%/RHc 90%）下的性能输出。而在干工况（RHa 50%/RHc 50%）下，工作温度对燃料电池输出性能的影响有限。

（2） 工作温度从70℃提高至90℃时，湿度变化引起的反应气体浓度变化越明显，各个温度获得最佳功率密度所对应的RH工况不同。相比于90℃工况下可在高阳极/低阴极工况（RHa 90%/RHc 50%）下获得最高的电流密度，燃料电池倾向于在70℃且高湿度的环境下获得更好的输出性能。

（3） 基于粒子群算法，优化得出当阳极RH为90%且阴极RH为0~50%范围内时，温度越高，获得最大功率密度的阴极湿度越小。工作温度和阴极RH分别控制在90℃和14.6%时可获得最高的功率密度，约为0.88 W/cm²。