质子交换膜燃料电池建模研究评述

摘要

车用动态工况下，质子交换膜燃料电池（PEMFC）内部压力、温度和气体浓度直接影响整个系统的性能和耐久性，因此需要借助数学模型分析电池内部传质、传热和化学反应过程，以优化燃料电池系统设计及控制。然而，燃料电池是一个非线性、多输入、强耦合的系统，模型的建立涉及从原子分子级尺度，到材料晶格结构级尺度，再到部件、电堆和系统级尺度。针对多尺度建模问题，从关键部件、单体、电堆及系统等角度，对不同尺度的建模方法进行全面梳理，并阐述了未来建模的挑战及发展趋势。

关键词

质子交换膜燃料电池; 多尺度; 建模

质子交换膜燃料电池（PEMFC）因其效率高、工作温度低和零排放等优点而成为新能源汽车的理想动力源，但其商业化应用仍然受到耐久性、可靠性等问题的限制。现有研究表明，车辆动态工况是导致车载燃料电池寿命衰退的主要原因^[1]。在动态工况下，燃料电池多维度的内部状态，如流道、气体扩散层内的温度、压力、气体浓度等波动较大，不当的操作条件会导致电池内部发生全局或局部水淹、饥饿现象。因此，为了提高燃料电池及其系统性能与寿命，需要建立全面且经过充分验证的数学模型，揭示燃料电池内部流体流动、传热和化学反应的多物理量多重耦合规律，为燃料电池及其系统的结构、材料优化和控制提供指导。

燃料电池内部过程涉及水、热、电多域问题，针对燃料电池建模，应当考虑燃料电池各反应、传质之间的相互耦合作用^[2]。此外，燃料电池内部过程涉及到不同的空间尺度和时间尺度^[3]。例如，电化学反应过程由催化剂分子和电解质的化学物理相互作用决定，在各个原子/分子的空间和时间尺度上呈现；电池内晶界传质倾向于发生在晶粒结构的空间和时间尺度上，涉及材料、界面、形态的化学和物理变化；而质量、电荷和热传递等传质过程则涉及从材料晶格到系统级尺度。不同的尺度特性，采用的建模方法也不相同。针对电子、分子等微观尺度，通常采用密度泛函理论等第一性原理；连续体、单体以及系统级等宏观尺度通常采用多个守恒方程，利用有限体积、有限差分等方法进行求解；微观和宏观之间的介观尺度建模可采用格子波尔兹曼等离散形态模型。据作者了解，目前没有一种模型适用于所有特征尺度，且每个特征尺度模型的计算时间也不相同。因此对PEMFC建模除了准确描述所有相关的物理化学机制之外，另一挑战是如何将这些机制在很大范围的空间和时间尺度上进行适当耦合^[4]。只有基于不同的尺度特性，才能更好地理解燃料电池内部详细的反应传质过程。因此，本文围绕质子交换膜燃料电池建模方法展开较为详细的评述，旨在梳理总结已有的建模方法，提炼出建模关键共性问题，促进先进燃料电池模型的开发与工程化运用。

1 燃料电池部件建模

PEMFC通常由质子交换膜（proton exchange membrane, PEM）、催化剂层（catalyst layer, CL）、气体扩散层（gas diffusion layer, GDL）和双极板（Bipolar plate, BP）等部件组成，为建立详细的燃料电池多域多尺度模型，需要从这些关键部件入手，深入了解各部件的相关机制，研究每个部件不同复杂程度的建模方法。

1.1　质子交换膜

车用燃料电池最常用的交换膜为全氟磺酸聚合物交换膜，例如Nafion系列交换膜，其导电性与膜内含水量密切相关，因此有必要针对交换膜内水传质特性进行研究。

（1）宏观传质建模

交换膜内含水量可由量纲为一符号$\lambda$表示，其定义为聚合物中每个磺酸分子所吸附的水分子数。理论上，在相同的温度和压力条件下，膜一侧与饱和水蒸气接触或与液相接触时，膜内含水量应相同。但有实验表明，虽然上述两种条件下的膜边界水活度相同，但是液相平衡水含量要高的多，该现象称为施罗德驳论^[5]。文献[6]采用热力学模型对其进行解释并对平衡水合度建模，其研究表明，膜在饱和气相中吸水时，首先是气相水在膜微孔内冷凝，该过程需要克服水的表面张力及膜内微孔的毛细压力，然后才发生液相水合作用。为了模型简便，大部分平衡含水量建模采用文献[7]提出的经验公式。近期，为了提高数值计算稳定性，文献[8]采用连续表达式来代替文献[7]中的分段式。

平衡吸附模型通常假设膜和相邻流体之间存在热力学平衡，然而燃料电池在变负载条件下运行时，电池内的气体相对湿度是快速变化的，因此需要描述膜内水传质的动态特性。水在膜内通常被认为是“溶解”相，主要存在3种传质机制：电子拖拽、浓差扩散和对流扩散。

第1种建模方法，假设膜内水传质由电子拖拽和浓差扩散主导，其守恒控制方程如下^[9]

$$\frac{{\rho }_{\mathrm{m}}}{E}\frac{\partial }{\partial t}\left(\omega \lambda \right)+\nabla \cdot \left(n_{\mathrm{d}}\frac{J}{F}\right)=\frac{{\rho }_{\mathrm{m}}}{E}\nabla \cdot \left(D_{\mathrm{m},\mathrm{d}}\nabla \lambda \right)$$

(1)

式中：$\omega$为聚合物体积分数；$E$为交换膜等效质量；$n_{\mathrm{d}}$为电子拖拽系数；$F$为法拉第常数；$J$为电流密度矢量；$D_{\mathrm{m},\mathrm{d}}$为膜内溶解水扩散系数；${\rho }_{\mathrm{m}}$为交换膜密度。

针对式(1)，目前多数文献采用平衡吸附假设，即膜吸水是瞬时发生，且水合作用后总是与催化剂层中的其他水相处于热力学平衡，此时膜内含水量的边界初始值由临近催化剂层内部的局部气相水活度和温度决定，且不包含源项。为保证催化剂层与膜界面处水浓度的连续性以及计算的简便性，可采用虚拟气相水浓度^[10]，将膜内溶解水相转换为气相，从而实现水传质的连续性，即式(1)可转换为

$$\frac{\partial }{\partial t}\left(\omega C_{\mathrm{m}}\right)+\nabla \cdot \left(n_{\mathrm{d}}\frac{J}{F}\right)=\nabla \cdot \left(D_{\mathrm{m}}\nabla C_{\mathrm{m}}\right)$$

(2)

式中：$C_{\mathrm{m}}$为膜内气相水浓度；$D_m$为膜内气相水扩散系数。然而，交换膜需要大约100~1000s才能达到气相水吸附平衡，20s达到80%程度的液相吸附平衡^[11]，因此瞬时平衡吸附假设不符合实际情况，故需对动态吸附进行建模。通常假设膜吸附/去吸附的水量是和局部含水量与平衡含水量之间的差值成比例^[12]，将其耦合到催化剂层内膜水传质方程的源项中，能描述动态吸附特性的影响。但该方面的研究仍然有许多特征未知，例如，吸附系数和去吸附系数在不同阶段的变化特性、电解质比表面积对吸附特性的影响以及膜内水去吸附后转变的是液相还是气相等问题，仍然是需要仔细考虑的。

第2种建模方法考虑对流效应为主导机制，氢离子传质采用Nernst-Planck方程描述，液态水动量由Schlögl方程计算^[13]。该模型假设膜为多孔介质，水主要由阴阳极压差驱动，因此在液态饱和条件下能够准确描述传质过程；但在较低含水量下，由于膜内离子簇不再连接良好，对流效应可忽略不计。另外，从物理角度来看，压力梯度驱动流过纳米级多孔膜的概念并不合理。

第3种建模方法不考虑驱动模式，将交换膜处理成尘埃粒子并引入扩散混合物中，使用广义Stefan‐Maxwell模拟混合物传质过程^[14]，Schlögl方程计算流体速度。但该建模方法假设较多，且模型参数难以获取，一般使用较少。故式(1)和式(2)仍然是目前交换膜内宏观水传质建模的常用方程。

（2）微观传质建模

为满足车用燃料电池的长寿命需求，PEM材料必须具有良好的抗氧化性、热稳定性以及机械稳定性。因此，为了进一步优化材料，需要从微观尺度进行建模。

分子动力学（Molecular dynamics, MD）能够根据预先定义的电场求解系统中所有粒子的牛顿运动，广泛运用于交换膜内微观传质分析。基于经典MD模拟可以发现当交换膜发生水合作用时，被吸收的水分子会与磺酸根结合，当膜吸收更多的水分时，水和酸性位点之间的吸引力较弱，位于磺酸位点附近的水相相对灵活，容易形成质子传递通道，从而增加了质子传导率^[15]。除了含水量外，Nafion侧链的长度以及等效质量对膜内质子、溶解相传质的影响仍可通过MD分析^[16]。质子交换膜内质子传递机理主要包括运载机理和Grotthuss机理，Grotthuss机理涉及水合氢和水的氢键重新排列，因此对该机制的理论研究需要具有捕获键断裂和形成的能力。经典MD无法处理键断合变化，因此无法描述Grotthuss机理。从头算分子动力学（Ab initio molecular dynamics, AIMD）能够描述该机制，但针对与质子扩散相关的长度尺度，计算时间过长^[17]。因此，文献[18]采用多态价键方法和反应分子动力学来准确描述PEM中的质子运载机理和Grotthuss机理, 与计算上昂贵的AIMD模拟相比，该组合模型可以在更长的时间和长度尺度上研究质子扩散。

为进一步研究非平衡条件下的质子传递，长度和时间尺度需要扩展到介观尺度。具有粗粒度模式的耗散粒子动力学（Dissipative particle dynamics, DPD）模型能够模拟数十纳米尺度下亲水憎水相的形态，该模型将水合离聚物进行分解，重组成介观尺度下的粗粒子，基于Flory-Huggins参数描述各种DPD颗粒之间的相互作用。然而由于自由度的减少，粗粒度模型无法准确描述与时间密切相关的物理特性，因此可能无法准确描述质子的动态传递。近期，文献[19]通过平滑粒子流体动力学方法，将系统视为准连续体，利用Nernst-Planck连续性方程模拟非平衡条件下质子和水的传质过程。该方法的优点在于系统被离散化为准粒子，每个粒子代表空间中的有限体积，通过求解密度场的时间演化或每个粒子的特定性质，可使用类似于经典MD的算法求解系统的连续性方程。

1.2　催化剂层

催化剂层建模主要考虑多孔电极内反应物和生成物的传质过程以及发生的电化学反应。

（1）传质建模

描述催化剂层最简单的方法是将其看成为膜与GDL之间的薄界面，但该方法无法描述多孔介质内部的传质现象。而宏观均质模型是目前催化剂层传质建模应用最广泛的方法，该类模型假设催化剂层内部各向同性，物理结构通过孔隙率、接触角和渗透率等参数度量。在多孔电极中，扩散占据主导机制，多相物质传质可基于修正的扩散系数利用Fick定律描述

$$\frac{\partial }{\partial t}\left(\epsilon C_i\right)+\nabla \cdot \left(u_{\mathrm{g}}{C}_i\right)=\nabla \cdot \left(D_{\mathrm{e}\mathrm{f}\mathrm{f},i}\nabla C_i\right)+S_i$$

(3)

式中：$C_i$为各气体浓度；$i$为组分代码，表示阴阳极的反应气体；$u_{\mathrm{g}}$为混合气相速度；$\epsilon$为孔隙率； $D_{\mathrm{e}\mathrm{f}\mathrm{f},i}$为有效扩散率；$S_i$为源项，包含电化学反应、膜内水的吸收和去吸收以及冷凝蒸发过程。

通常使用合适的自扩散或二元扩散系数的Fick定律是合理的，但为了更加严谨，有作者采用Stefan-Maxwell方程来描述扩散过程^[20]。气体速度和压力主要基于动量守恒方程和质量守恒方程联合求解。在两相流模型中，表面张力在多孔电极中占主导地位，动量守恒方程可简化成两相Darcy方程。与单相模型方程相比，两相模型采用了相对渗透率，表示由于存在另一相而减小的孔隙空间对每个相的影响。然而实际上，CL的内部孔隙率、渗透率和相对电导率等特性存在各向异性，且宏观均质模型往往忽略了由膜两端电荷聚集产生的双电层电容效应，因此宏观均质模型无法详细描述CL的局部传质^[4]。

目前团聚体模型能够代表真实的CL结构而得到广泛认可，通常认为存在三类团聚体，矩形团聚体、圆柱形团聚体和球形团聚体。根据扫描电子显微镜结果，球形团聚体模型是反映真实CL结构最合适的模型^[21]。在球形团聚体中，假设CL由大量微米级的球形聚合物连接组成，其中填充铂催化剂颗粒和电解质簇，因此CL的孔隙率和电解质分数可通过铂比率和团聚体半径等参数计算，而非直接选取经验常值。文献[22]对经典的团聚体模型进行了改进，团聚体内部催化剂颗粒采用的是离散随机分布，而不是理想的均匀混合分布。

为了分析团聚体之间孔隙形成的“次级”孔内的传质过程，可采用直接数值模拟（Direct numerical simulation, DNS）模型。基于电子显微镜图像，可将CL理想化为具有一定规律或者随机分布的微孔结构，直接利用守恒方程求解^[23]。然而DNS模型需要精确测量CL的物理结构参数，且对CL结构进行重构会引起额外误差；此外，DNS模型通常是单相模型，无法考虑液态水分布问题。孔隙网络模型（Pore Network Modeling, PNM）常用于多孔介质内的两相传质计算，并能够将局部特征与其有效传输特性联系起来。在PNM中，次级孔隙理想化为由涂有薄Nafion层圆柱形喉部连接的规则网络^[24]，喉孔通过没有体积的节点相互连接。由于CL内部环境比GDL更加复杂，PNM模型通常忽略CL中团聚体内和电解质内的传质，并不与电化学反应模型耦合。因此，文献[25]提出了一种新的PNM模型，该模型耦合了质子、氧气传递以及电化学反应方程，能够分析CL的结构参数对燃料电池电压特性的影响。

格子玻尔兹曼方法（Lattice Boltzmann method, LBM）同样可分析多孔介质内的传质过程。与DPD模型中的随机碰撞和MD模型中的单个分子追踪不同，LBM通过粒子在规则晶格上移动和碰撞来模拟流体传动，每个粒子由高度灵活和数值稳定的多个分子组成^[26]。然而在CL中，LBM与电化学方程难以耦合，因此为了实现耦合计算，CL被视为薄界面^[27]或只有固相而没有电解质相^[28]。近期，文献[29]首次将LBM与团聚体模型结合，氧气传质以及氮气传质通过LBM团聚体模型解释。

（2）电化学反应建模

阴极的氧化还原反应（Oxygen reduction reaction, ORR）由于过电势而限制了燃料电池的功率输出，因此需要在原子尺度分析ORR的反应机制。ORR是通过2电子转移生成过氧化氢中间产物再生成水，还是通过4电子转移直接生成水，仍然存在争议。密度泛函理论（Density functional theory, DFT）可用于模拟吸附材料和基材之间的相互作用，因此广泛用于ORR反应机理研究。基于DFT模拟可发现，由于燃料电池处于酸性环境，O-O键不先发生断裂，而是与H生成OOH基团^[30]。但反应步骤机制受催化剂层表面OH覆盖率的影响，有作者表明当覆盖率小于1/3时，优先发生2电子转移反应；反之，则优先发生4电子转移反应^[31]。溶剂和催化剂表面吸附物质会影响反应物分子的吸附，因此还需分析电化学反应中原子尺度的吸附动力学。文献[32]基于MD模型分析了O₂和H₂O分子在3种不同催化剂表面的吸附方式和吸附能，并还研究了温度对ORR反应的影响。蒙特卡罗动力学（Kinetic Monte Carlo, KMC）同样可描述电化学反应涉及的表面动力学、覆盖物演化以及吸附、解吸和反应速率等机制。与DFT第一性原理相比，KMC能够描述更宽的尺度范围，但KMC模型依赖催化剂表面上的反应物和中间产物的自由能计算。因此，文献[33]结合KMC和DFT方法，KMC建模所需的参数由DFT方法提供。近期，有作者将KMC模型与连续传质模型进行耦合，可模拟介观尺度下的催化剂表面化学反应过程^[34]。

燃料电池电化学反应具有多个复杂的步骤,为了计算简便，通常采用Butler-Volmer方程描述阴阳极氧化还原反应速率，并与质子、电子守恒方程的源项进行耦合，共同描述燃料电池内质子、电子传质过程^[21]。在宏观均质CL模型中，简化的反应气体扩散会低估高电流密度下的浓差损耗。然而在团聚体模型中，氧气首先通过催化剂层的孔隙扩散到电解质膜表面，然后溶解到电解质中，孔和电解质膜界面处的氧气浓度通过Henry定律描述。氧气在电解质膜内传质期间，由于没有发生电化学反应，因此可通过Fick定律描述溶解通量。最后，氧气到达催化剂表面并发生电化学反应，更详细的传质过程可参考文献[22]。由于团聚体模型描述了更详细的氧气传质和溶解过程，故可捕获纳米尺度下颗粒相互作用引起的扩散损失，而且与宏观传质模型耦合方便，计算量小。因此，建议宏观尺度模型的电化学反应建模采用团聚体模型，但与团聚体模型相关的材料参数需要仔细考虑。

1.3　气体扩散层

气体扩散层为催化剂层提供结构支撑，并对反应物、产物以及离子进行传输，对其建模主要研究其内部两相多组分流动、热量和质量传递耦合问题。

基于宏观连续假设的GDL内部两相传质建模主要有两种方法：多相混合模型（Multiphase mixture model, MM）^[35]和双流体模型^[36]。MM模型将气相水和液相水在质量和动量守恒方程中作为混合物进行求解，然后基于混合物平均特性获取液态水分布。在MM模型中，由于气液相共用守恒方程，模型复杂度和计算量降低了，但通常假设气相和液相之间保持热力学平衡，因此动量和组分守恒方程的源项中不包含水相变过程。与MM模型相比，双流体模型中添加了额外的液相水传输控制方程，并可考虑水的相变速率对液态水传质的影响。在GDL多孔介质中，孔壁的黏附效应变得显著，多孔毛细管效应和表面润湿性对液态水传质具有较大影响，而诸如重力的体积力不会显著影响其动量，因此，基于广义达西定律和毛细管压力函数，多孔介质中液体传质的守恒方程可表示为^[24]

$$\begin{array}{l}\frac{\partial }{\partial t}\left({\rho }_{\mathrm{l}}\epsilon s\right)+\nabla \left(\frac{K_{\mathrm{l}}{\mu }_{\mathrm{g}}}{K_{\mathrm{g}}{\mu }_{\mathrm{l}}}{\rho }_{\mathrm{l}}{u}_{\mathrm{g}}\right)=\\ -\nabla \left({\rho }_{\mathrm{l}}\frac{K{K}_{\mathrm{l}}}{{\mu }_{\mathrm{l}}}\frac{\mathrm{d}{\mathrm{P}}_{\mathrm{c}}}{\mathrm{d}\mathrm{s}}\nabla s\right)+S_{\mathrm{l}}\end{array}$$

(4)

式中：${\mu }_{\mathrm{g}}$和${\mu }_{\mathrm{l}}$分别为气液相黏度；$s$为饱和度；$K$为相对渗透率；$K_{\mathrm{l}}$和$K_{\mathrm{g}}$分别为液体渗透率和气体渗透率；$S_{\mathrm{l}}$代表电化学反应、蒸发和冷凝过程中的质量源项；${\rho }_{\mathrm{l}}$为液态水密度；$P_{\mathrm{c}}$为渗透压力。从上式可知，液相饱和度的计算与毛细管压力函数以及相对渗透率大小密切相关。目前，Leverett方法广泛用于PEM燃料电池液相水传质的毛细管渗透压力建模^[37]。但GDL内部传质受到多孔介质湿润性分布的影响，针对高度异性的薄多孔材料，Leverett方法精度不够。因此，需要考虑温度、压缩性和润湿性，构建新的渗透压力复合函数，具体函数研究现状可参考文献[37]。另外，为了描述冷凝和蒸发相变过程，在宏观尺度上，可由饱和蒸气压力和水蒸气分压的压差描述。

虽然上述宏观模型可以描述GDL内两相流动，但由于孔径与GDL厚度之间缺乏显著差异，因此宏观模型无法描述局部传质特性，需借助如前所述的PNM和LBM介观尺度模型。在PNM中，GDL同样被理想化为遵循一定尺寸分布的孔体和喉孔的规则立方网络，通过考虑孔隙中的毛细管力和通过喉孔的入口压力、黏滞压力来获取GDL中两相传质过程和液态水分布，该模型相对宏观连续模型来说具有更高的精度。然而水进入喉孔条件对液体通过GDL的行为具有重要影响，故应研究更具代表性的进入边界条件，然而该边界条件与CL中水的产生、传输以及通过微孔层的过程密切相关，使得边界条件难以确定，这是PNM模型使用的主要难点^[38]。相比PNM模型，LBM模型具有更高的灵活性，能够模拟GDL内部复杂结构的流体流动，渗透率参数可通过该方法计算得出^[39]。由于PNM假设孔隙仅由一个网格点表示，相比之下，LBM假设每个空间方向上的每个孔径使用多个网格点描述，因此LBM具有更高的精度，但计算时间和资源需求也更多。

1.4　双极板

双极板内的气体流道（gas channel, GC）影响反应物和生成物的传输以及反应物的利用率，对其建模主要研究不同操作条件和流道物理特性对流道内部及GDL与流道界面的两相流动和液滴动态特性的影响。

气体流道内传质主要为对流效应，通常假设流体为理想流体，采用Navier-Stokes方程描述其传质过程。但液滴或水膜在气体流道中会发生局部集中、分散和累积现象，这些模型无法识别液态水的动态行为；而PNM和LBM不适用于模拟燃料电池流道内宏观层面的两相流动。为此，需要采用界面跟踪的建模方法。流体体积法（Volume of Fluid, VOF）能够通过求解单组动量方程来求解整个流道内每种流体的体积分数并捕获液滴形状。此外，该模型还能考虑液滴表面张力和壁面黏附效应，因此广泛用于分析气体扩散层表面湿润性、通道壁面湿润性、空气入口流速和注水位置、注水速度等条件对气体流道内水滴动态特性的影响^[40]。

通常VOF方法的计算步长时间较短，约为10^-6~10^-3s，而流道内形成段塞流可能需要几分钟甚至几小时，这使得VOF方法也很难与其他的传输现象耦合起来，因此大部分VOF方法只应用于流道内，即其他的传质现象被忽略。因此文献[41]建立了与VOF方法耦合的三维动态多相燃料电池模型，其中VOF运用于通道和多孔介质内，并将流体流动、热量传导、电子和质子传输与电化学反应相耦合，分析了两相流动模式以及液态水动态特性对速度场、压力分布、反应物浓度、电流密度和温度分布的影响。但耦合模型为了缩短计算时间，对液滴在流道的初始位置和注入速率进行了简化，因此，只适用于特定的条件。目前，一般建立燃料电池单体、电堆尺度数值模型时，往往忽略流道内液滴的动态特性，或者直接将流道内的水当作单相处理。

2 燃料电池单体、电堆及其系统建模

研究各部件不同尺度建模方法后，还需将前面所述各部件的建模方法进行耦合，从而建立完整的燃料电池数值模型，以分析各部件特性对整个单池、电堆性能的影响。此外，由于燃料电池的控制策略必须考虑重要执行器和外围系统的动态特性，因此有必要研究系统级的建模问题。

2.1　燃料电池单体、电堆

燃料电池单体、电堆的数学模型可分为分布参数模型和集总参数模型。

（1）分布参数模型

如前所述，燃料电池各个部件都有对应尺度的建模方法，目前常用的燃料电池单体、电堆分布参数建模仍主要在宏观尺度层面，即将所有多孔电极处理成同性的多孔介质或简单的异性，通过耦合质量、动量、组分和电荷守恒等方程，计算电池内部压力、速度、组分浓度和液态水等空间分布^[3,4]。但是，多数文献假设燃料电池内部温度恒定，然而温度对燃料电池内部质量、电荷传递和相变过程具有较大影响，从而影响电池输出性能；反之内部传质、相变过程也会影响温度变化。为此，有必要建立燃料电池热模型 ^[42]。此外，目前常用的计算流体动力学（Computational Fluid Dynamics, CFD）方法有有限差分法、有限元法和有限体积法^[12]，其中有限体积法运用最广泛，目前多种CFD商业化软件封装了上述计算方法。但燃料电池分布参数模型具有非线性强、耦合性强、求解计算量大等特点，因此有必要研究更高效的数值计算方法。针对模型计算域问题，根据计算维度通常可分为一维、二维、三维模型。目前对于三维燃料电池单体分布参数模型，通常是在单流道水平面上进行模拟，虽然计算域包含了PEM、CL、GDL和GC等部件，但大部分模型没有包含冷却流道内的仿真。然而燃料电池单体温度梯度分布对水热管理的影响不可忽略，因此在三维模拟仿真也应考虑冷却流道的影响^[43]。

宏观均质模型在计算上是高效的，并且可以处理大规模的模拟问题，但无法描述多孔介质内细节的局部运输，而前面提到的MD、KMC、LBM、PNM、团聚体模型等能够有效描述燃料电池局部特性，但除了团聚体模型外，其他方法需要大量计算资源，且计算域过小，此外不同尺度的计算方法难以直接进行耦合计算。为此需要采用间接耦合方法，在催化剂层面，可利用KMC模型计算催化剂上的瞬时反应速率，反应物浓度可利用CFD传质模型求解，KMC模型所需的反应物及中间产物活化能参数可由DFT计算^[34]，这种通过信息交互的间接耦合方法使得多尺度建模成为可能。此外，文献[26,27,28] 也研究了有限体积法和LBM，有限体积法和MD，有限体积法和KMC的耦合计算。但耦合模型在不同计算域界面之间信息交换困难，而且依托了大量假设和简化，实际应用会受到一定限制。

（2）集总参数模型

集总参数模型不考虑燃料电池内部详细机制，忽略空间上的变化，基于一定的经验公式，在水平方向上分析其动态特性，因此该类模型通常称为零维模型，具有仿真用时短、操作方便等优点。文献[4,42]总结了目前常用的控制导向动态集总参数模型。但是，集总参数模型越复杂，涉及到的未知参数越多，还需利用实验得到的离线数据，对未知参数进行寻优，燃料电池模型参数识别方法可参考文献[44]。另外，集总参数模型相关变量都是“平均值”，例如，代入电压方程的阴阳极气体分压都是平均压力。然而燃料电池内部状态具有非均匀性，因此需要对分布式参数模型进行降维、降阶处理，开发能够描述参数分布信息的集总参数模型。文献[45]对质子交换膜燃料电池单体沿质子传递方向上各个部件内的多状态水传递过程进行建模，通过数学解析方法推导了两相水一维稳态分布的解析解，该简化模型计算简便又能够描述空间分布信息。具有分布式特性的集总参数模型是以后控制导向模型的发展趋势。

2.2　燃料电池系统

燃料电池系统建模在系统设计和开发以及燃料电池控制策略的优化和测试中起到决定性作用。由于实际系统的控制策略必须考虑重要执行器和外围系统动态特性，因此必须建立适用于控制的电压、供气系统和热管理系统等子系统的动态模型。文献[46]建立了一个包含空压机、歧管流动、流道的气体动态传输及膜干湿变换的集总参数动态模型，该模型已得到广泛认可，目前许多用于气体供给系统控制的模型均以此为参考，但该模型忽略了电堆的温度变化。文献[47]除了对燃料电池空气供应子系统建模外，还集成了氢循环的引射器模型。车用燃料电池系统还需通过DC/DC耦合连接于高压母线上，为此也需要集成燃料电池系统模型与DC/DC模型，以分析DC/DC开关频率、控制模式对燃料电池性能和内部状态的影响^[48]。为进一步探究辅助单元动态特性对电堆内部状态分布的影响，文献[49]建立协同闭环仿真模型，利用集总参数模型描述辅助系统的动态特性，并为电堆数值模型提供动态边界条件；电堆分布参数模型分析电堆内部状态的动态分布，并为辅助系统模型提供相应的变量反馈，该方法为以后的多尺度仿真提供了新的思路。神经网络、支持向量机等智能算法也可用于燃料电池系统建模^[50]。但智能算法模型通常只研究系统输入和输出之间的响应特性，而不关心系统本身内部的耦合特性，所以不能反应电池系统内部状态。但在未来大数据趋势下，结合数据驱动和机理模型也是未来系统级模型的方向，可用于燃料电池故障诊断研究。

2.3　模型验证

为说明模型仿真结果的正确性，需要对模型进行验证。现阶段，几乎所有文献都采用极化曲线对其外特性验证。然而极化曲线受到多个参数影响，两组不同的参数可能会导致几乎相同的极化曲线。为了进一步验证模型合理性，建议在相同的极化曲线中将各种电压损失与实验数据进行对比；同时，应将两个或多个不同参数或不同操作条件下的极化曲线仿真结果与实验数据进行对比。燃料电池交流阻抗谱能够从频域的角度获取电池内部状态信息，利用快速傅里叶变换或小波变化，将模型输出的时域特性换成频域特征，也能够对模型进行进一步验证，相比传统的极化曲线更有说服力，而且阻抗谱获取也较容易。

除了利用外特性数据，另一方面可直接利用电池内部测量结果对模型仿真结果进行验证。燃料电池内部测量实验侧重于内部液态水和冰的分布测量，例如核磁共振技术、中子成像技术、扫描电镜和X射线技术，但这些设备昂贵，操作复杂，实验技术要求也较高。为此，高校更倾向于采用技术成本较低的光学成像法；将燃料电池双极板改成透明材料，或在断面部分开窗，以获取内部可视通道，利用光学照相设备直接观察内部液态水分布^[51]。

除了利用实测液态水分布对数值模型验证外，还需测量电池内部气体分布，可分为直接测量和间接测量。直接测量法需对双极板流道区域开孔设计，方便将毛细管插入孔内，对流道内气体进行采样，利用质谱仪分析各气体分压，但进入质谱仪的采样气体需加热至100℃以上^[52]。间接测量主要是对电流密度分布进行测量，根据电池输入信息反推评估出各气体浓度分布，目前主要分为电阻网络技术、印刷电路板技术和电磁测试技术^[53]。除了电流密度分布测试外，还可利用分区技术实现分区阻抗测量，进一步反映内部状态分布信息^[54]。温度分布测量一般通过在燃料电池不同位置植入温度传感器实现，传感器类型有热电偶、热敏电阻和光纤光栅等，其他测量方法包括红外热成像法和激光吸收光谱法等。

3 质子交换膜燃料电池建模挑战

基于上述质子交换膜燃料电池建模方法评述，目前建模领域仍然面临如下挑战：

（1）精确两相数值建模

针对PEMFC内两相流建模，传统的混合模型和双流体模型采用宏观均质假设，依赖材料和传质过程的经验参数，这种均质近似的准确性需要进一步研究；能够捕获液滴动态行为的宏观VOF方法由于计算时间过长，无法完全考虑电极的微观结构；介观建模方法（例如LBM）可对多孔介质中的两相特性建模，但难以扩展至大尺度仿真（例如PEMFC单体、电堆）中，且难以与电化学反应耦合。因此，建立精确的多尺度两相流模型，非常具有挑战。

据此，为开发更好的两相流模型，需要为CL和GDL定义易于测量且严格有意义的传输参数；此外，需要结合原位和非原位等实验手段以更准确地了解燃料电池中实际发生的两相流动现象（流动模式、压降、流体分布、流动迟滞）并验证当前的数值模型；同时，模型研究应该结合宏观、介观和微观尺度模型的优点，同时避免各自的缺点，以模拟燃料电池中的多尺度传质现象。另外，需开发更高效的数值计算方法，以缩短宏观尺度和微观尺度的计算能力差异。

（2）精确的热学建模

针对热学建模，通常会采用如下假设：低电流密度下或小型实验电池采用绝热假设；热接触内阻被忽略；仅考虑气流的对流传热；气流与双极板之间不存在热交换；忽略冷却液流道相互作用；双极板不会与环境进行热交换。这些假设会限制温度梯度分布的准确性，关于不同电池部件的热导率也缺乏准确的实验数据；此外，多数文章分析热传递过程时没考虑两相流传质，反之亦然。

据此，需要了解电池各部件准确的热导率，特别是从各部件的多维和各向异性等角度获取；此外，需充分考虑上述假设合理性，特别是针对大面积车用燃料电池，其产热功率较大且热交换过程不可忽略；同时，燃料电池的热学分析需结合两相流、质量和电荷传质过程。

（3）冷启动数值模型不确定性

冷启动数值模型对分析燃料电池内部结冰融冰和冷启动策略具有重要意义，但模型仍存在着许多不确定性。多数文献假设膜结合水是从过饱和状态直接瞬时结冰，而这种假设与实际相变过程存在差异；此外，过冷水通常假设为膜结合水和冰之间的中间产物，然而过冷水是具有随机性的，它的存在目前也仍存在争议。在冷启动模型验证中，常用的变量包括电压输出、累积水量和电池温度演变。而微观实验技术（如扫描电镜）虽然能够显示催化剂层内冰的纳米级分布，但由于其静态性质，该方法无法捕获液体到冰的转变过程；其他一些技术（例如中子成像）可以捕获冰的形成动力学，但是无法在纳米级测量中进行观察；现有实验技术导致我们对冷启动结冰机理的认知不充分，致使模型仿真结果存在不确定性。

据此，未来应探索更先进的观测技术，以便深入了解结冰机理；此外，与两相流模型类似，需要考虑多尺度建模特性，尤其多孔介质内冰的形成特征；同时，通过开展冷启动性能的系统级实验以及可以验证材料特性和基础理论的组件级实验来进行模型验证，降低后续建模过程中的不确定性。

（4）多域多尺度工程级电堆模型

如前所述，燃料电池内部涉及气、水、电、热等多域问题，但对工程级电堆来说，应力作用对电堆输出性能的影响同样不可忽略。系统装配、温度变化、变载工况和结冰融冰等会导致电堆内部应力发生变化，从而影响电堆内部特性，然而目前针对应力场耦合建模问题研究较少。此外，工程应用电堆的单体数量多至上百片，仅仅研究一片单体或几片是不够的，但片数过多则会导致计算需求过大，与现有计算能力和资源矛盾。

据此，为建立工程应用的电堆模型，需利用应力对于材料结构参数的影响构造应力与动量、组分、电荷传输和电化学过程以及热传输过程的弱耦合关系，迭代求解。此外，对工程电堆内部物理过程进行简化，提取关键过程因素，简化边缘物理过程；分析网格和计算精度之间的敏感性，从而达到削减网格，减少计算资源需求的目的。

（5）燃料电池系统模型实用性和完整性

虽然，目前关于系统级模型展开了大量研究，但是模型的实车运用能力，仍然需要进一步优化。例如，文献[46]中空压机的控制量为驱动电压，而实际的控制量为转速，且模型中空压机转动惯量、压缩效率等参数不易获取，所以目前车用控制器几乎都采用查表的方式进行控制。另外，目前系统级模型缺乏完整性，例如空气侧建模通常不考虑加湿器，而车用加湿器通常采用膜加湿器，利用电堆阴极排气对进气加湿，这种排气和进气之间的加湿耦合关系，目前相关建模研究较少；此外，氢气循环回路建模缺乏研究。

据此，应当对燃料电池系统重要的执行器进行测试，并建立参数容易获取、具有一定物理意义的模型，使其可应用于实车控制器上；另外需要从燃料电池系统的子系统出发，忽略弱耦合的状态变量和控制变量，建立完整的燃料电池系统模型。

4 结语

本文评述了质子交换膜燃料电池建模研究现状。在燃料电池部件方面，重点介绍了交换膜、催化剂层、气体扩散层、气体流道的宏观至微观等不同尺度建模方法；在燃料电池单体、电堆和系统层面，总结了目前常用的模型类型和适用场景；在燃料电池模型验证方面，总结了目前常用的模型验证方法以及实验测量手段。最后对现有建模研究的挑战和未来可能解决途径进行了阐述。虽然取得大量的理论和实验成果，但质子交换膜燃料电池建模研究还处于不断发展和完善阶段，未来仍然需要在以下方面进行拓展：

（1）结合宏观和微观建模方法，构建更精确的两相模型来模拟从CL表面上水的形成到GDL和气体流道中的冷凝和累积以及从气体流道中去除的迁移过程，以优化燃料电池内部水管理能力。

（2）分析应力对燃料电池内部特性影响，构建包含“水-热-电-气-力”等多域的多尺度模型，以充分揭示燃料电池内部多物理量多参数耦合关系。

（3）优化工程应用电堆的数值模型建模方法，简化边缘物理过程，分析网格和计算精度之间的敏感性，削减网格数量和计算资源需求，实现工程电堆的数值建模。

（4）优化燃料电池系统模型，尤其是空压机、氢循环系统、膜加湿器等关键部件，使其能够满足现阶段市场上燃料电池系统控制器应用；此外，需开发适用控制器应用的电堆分布参数模型，为分布式内部状态观测器设计提供基础。

（5）探究先进的燃料电池可视化技术，从微观角度分析燃料电池内部相变瞬态过程。

（6）开发新的模型求解计算方法，缩短宏观、微观之间计算能力需求差异，减少计算资源需求。