韓麗輝,于春梅,馮根生

(北京科技大學 冶金與生態(tài)工程學院，北京 100083)

燃料電池是將燃料的化學能直接轉(zhuǎn)換成電能的裝置。燃料電池具有效率高、無污染、重量輕等優(yōu)點，被稱為21世紀的綠色環(huán)保能源，具有很大的發(fā)展?jié)摿蛻?yīng)用前景[1]。在燃料電池的研究中，實驗室的物理實驗和數(shù)值計算必不可少，可靠的計算程序不僅為物理實驗提供理論指導，還可以優(yōu)化設(shè)計燃料電池、節(jié)約實驗成本、提高科研工作效率。數(shù)值模擬軟件Fluent在解決流體計算問題方面應(yīng)用很廣泛，并且該軟件具有燃料電池模塊，但其質(zhì)子交換膜燃料電池模塊是針對氫氧燃料電池開發(fā)的，當用戶的燃料電池是甲醇或其他質(zhì)子交換膜燃料電池時，該模塊就不能完全適用，用戶必須使用UDF編寫自己的燃料電池程序。本文從Fluent數(shù)值求解機理、燃料電池問題基本設(shè)置、編寫用戶UDF、程序調(diào)試等方面，詳細闡述了如何利用Fluent軟件來數(shù)值求解質(zhì)子交換膜燃料電池問題。

1 Fluent求解機理

1.1 基本數(shù)據(jù)結(jié)構(gòu)

了解Fluent數(shù)據(jù)存取方式及程序求解機理是編寫UDF的前提和基礎(chǔ)。用戶自定義函數(shù)需要從Fluent求解器中讀取數(shù)據(jù)，數(shù)據(jù)的基本載體為網(wǎng)格。每個網(wǎng)格單元稱為cell，每個cell都有節(jié)點(node)、面(face)、網(wǎng)格中心點(cell center)，見圖1。每個網(wǎng)格都由面包圍而成。Thread 是一種數(shù)據(jù)結(jié)構(gòu)，可以用來存儲一組網(wǎng)格和一組面的信息。與網(wǎng)格相關(guān)的數(shù)據(jù)類型為Node、face_t、cell_t 、Thread、Domain，通過這些數(shù)據(jù)類型定位到數(shù)據(jù)載體。本文中求解燃料電池問題的UDF程序會用到Fluent中已定義的一些基本數(shù)據(jù)宏[2-4]見表1，如網(wǎng)格中心點坐標，網(wǎng)格體積、流體密度、流體溫度、流體壓力、混合物組分、變量梯度、多孔介質(zhì)孔隙度等。 表1中，ND ND為Fluent中全局變量，對于二維問題，該值為2，對于三維問題，該值為3；另外宏C_P(c,t)讀取到的壓力是表壓力(Gauge Pressure)，其絕對壓力等于該值加上操作壓力，F(xiàn)luent中全局操作壓力變量為op_pres。

圖1 Fluent 中基本數(shù)據(jù)結(jié)構(gòu)

表1 Fluent中部分基礎(chǔ)數(shù)據(jù)宏

1.2 方程求解過程

求解燃料電池問題需要求解自定義標量傳輸方程，F(xiàn)luent中以壓力為基本項的耦合求解器求解過程如圖2所示。在迭代計算外要執(zhí)行兩步運算：一是初始化，即用Fluent初始化界面中的數(shù)據(jù)來初始化求解變量，之后調(diào)用用戶自定義的初始化程序，后者數(shù)據(jù)會覆蓋前者數(shù)據(jù)；二是調(diào)用用戶自定義邊界條件程序(PROFILE UDF)。

求解迭代循環(huán)過程：首先從執(zhí)行用戶定義的調(diào)整程序(ADJUST UDF)開始，之后耦合求解連續(xù)方程和動量方程；最后順序求解能量方程、組分傳輸方程、湍流方程、用戶自定義標量傳輸方程；所有方程求解結(jié)束后執(zhí)行物料屬性更新程序(用戶自定義PROPERTY UDF)；最后檢查是否滿足收斂標準，不滿足進行下一次迭代計算，滿足則退出循環(huán)計算結(jié)束。

了解fluent控制方程求解過程對于用戶編寫自己的UDF非常重要，根據(jù)求解過程，用戶就能在合適的位置給自定義內(nèi)存賦值。用戶感興趣的計算數(shù)據(jù)是通過自定義內(nèi)存存儲的，自定義內(nèi)存的賦值位置可以放在程序調(diào)整宏 、控制方程源項、物料屬性UDF，而不能放在迭代循環(huán)外邊的初始化及邊界條件的UDF 程序中。

圖2 Fluent中壓力基本項求解過程

2 燃料電池問題的基本設(shè)置

2.1 調(diào)入氫氧燃料電池模塊

通過文本框命令行調(diào)用質(zhì)子交換膜燃料電池模塊。雖然Fluent提供的質(zhì)子交換膜燃料電池為氫氧燃料電池，不完全適用甲醇燃料電池問題，但同屬于質(zhì)子交換膜燃料電池，其共同之處是可以借用的，如燃料電池各層材料屬性、水合量擴散系數(shù)、水合量在催化層固定值等。

2.2 基本模型的選擇

自定義標量UDS及內(nèi)存UDM。在界面中定義自定義標量UDS及自定義內(nèi)存UDM個數(shù)，同時要保證界面中定義的數(shù)量大于UDF程序中需要的數(shù)量。在質(zhì)子交換膜燃料電池中至少需要3個標量，即電子電勢、離子電勢及水含量。自定義內(nèi)存?zhèn)€數(shù)根據(jù)用戶需要而定，至少包括體電流密度、電流密度矢量、水含量、電滲拖曳系數(shù)、反應(yīng)熱等。

選擇層流模型、能量方程及組分傳輸模型。在組分傳輸模型中不考慮體化學反應(yīng)，化學反應(yīng)是通過標量及源項的設(shè)置體現(xiàn)出來。

2.3 流體屬性的設(shè)置

從Fluent物料數(shù)據(jù)庫中選擇甲醇空氣燃燒混合物,該混合物包括甲醇、氧氣、水、二氧化碳、氮氣。按順序設(shè)定流體組分，其索引號從0開始，依次為1、2、3、4等，該索引號在UDF中定義物料組分擴散系數(shù)時會被使用。流體屬性除了密度、比熱、黏度等基本參數(shù)外，還需要設(shè)定2個參數(shù)，且這2個參數(shù)對計算結(jié)果有很大的影響。一個參數(shù)是流體組分擴散系數(shù)。在多孔介質(zhì)結(jié)構(gòu)中，氣體分子往往被毛細管壁所阻礙，為了計算這些阻礙的影響，修正擴散通量，一般是通過修正擴散系數(shù)實現(xiàn)的，該系數(shù)主要影響組分傳輸方程。在質(zhì)子交換膜燃料電池問題中計算組分擴散系數(shù)時可以使用稀釋近似模型[5]：

(1)

另一個參數(shù)是自定義標量擴散系數(shù)(UDS,diffusivity)。對于電子電勢和離子電勢來說，其擴散系數(shù)即在不同流體區(qū)域的電導率，交換膜的電導率與電池溫度及水含量關(guān)系[6]為

(2)

式中，σmem為交換膜的電導率(Ω-1·m-1)；β為膜的電導率系數(shù)；λ為膜中水含量；ω為膜的電導率指數(shù)。

擴散層、催化層的電導率為常數(shù)，集流板固相區(qū)域的電導率需要在固體材料屬性UDS diffusivity 中設(shè)定。

2.4 區(qū)域的設(shè)置

共有9個計算域，包括陰陽極集流板(2個)、陰陽極流體通道(2個)、陰陽極擴散層(2個)、陰陽極催化層(2個)、質(zhì)子交換膜(1個)。其中集流板為固體域，流體通道、擴散層、催化層、膜為流體域。區(qū)域?qū)傩缘脑O(shè)置包括多孔介質(zhì)、源項及固定值。膜、催化層和擴散層均為多孔介質(zhì)區(qū)域，要選擇各自的多孔介質(zhì)材料。多孔介質(zhì)主要引起動量損失，區(qū)域為多孔介質(zhì)域時，在動量方程中自動引入動量損失源項。計算該源項需要黏性損失系數(shù)、慣性損失系數(shù)及孔隙度。當流動是層流流動時，慣性損失可以忽略。燃料電池數(shù)值計算實際上是在計算流動問題基礎(chǔ)上加載電勢傳輸方程，同時由于電化學反應(yīng)引起能量方程、組分方程、連續(xù)方程的源項發(fā)生變化，必須要加載源項(UDF 程序)。流體穿過多孔介質(zhì)區(qū)域從陽極到陰極或從陰極到陽極，這種穿透運動有兩種驅(qū)動力，一是由于動量驅(qū)動的，二是由于濃度差驅(qū)動的，為了使兩極流體不相混合，消除動量驅(qū)動作用的影響，在膜域?qū)α黧w運動的速度設(shè)為固定值0。另外膜的水含量影響膜的導電率及水的電滲效果，所以水含量只存在于催化層和膜，對于其他流體域水含量都需設(shè)為固定值0。

2.5 邊界條件的設(shè)置

陰陽極入口的質(zhì)量流量、流體溫度，混合物中各組分質(zhì)量比例、出口壓力、出口回流氣體溫度、集流板溫度、操作壓力都是可變量，可根據(jù)情況設(shè)定。陰陽極集流板電勢或電勢梯度如下：

陽極電子電勢?an：

?an=0

(3)

陰極電子電勢?cat：

?cat=Vcell

(4)

膜離子電勢?m：

(5)

式中n為膜厚。

2.6 求解控制及初始值的設(shè)置

壓力與速度耦合求解采用SIMPLE算法，采用固定循環(huán)模式，即按變量順序求解，使用BGSTAB方法增加求解過程的穩(wěn)定性。給自定義標量電子電勢和離子電勢賦予合理的初始值，是計算是否順利進行的關(guān)鍵，如果電勢初始值不合理，自定義標量殘差曲線出現(xiàn)振蕩或計算出現(xiàn)錯誤而停止。合理的電勢初始值能夠使得陰陽極都有合理的過電勢，即電子和離子運動的驅(qū)動力，這樣才能使電化學反應(yīng)連續(xù)發(fā)生。一般離子電勢的初始值為負數(shù)，如-0.1。

3 編寫燃料電池程序

在Fluent軟件提供的動量、能量、連續(xù)、組分傳輸?shù)确匠痰幕A(chǔ)上加載電子電勢及離子電勢標量傳輸方程，并對原有方程加載源項，并修改混合流體的相關(guān)屬性，最后數(shù)值求解方程組，得到流場、質(zhì)量組分、電流密度、溫度場的分布情況。

3.1 質(zhì)子交換膜燃料電池基礎(chǔ)

質(zhì)子交換膜燃料電池(PEMFC)結(jié)構(gòu)如圖3所示，其陰陽兩極各包括集流板(collector)、流道(gas channal)、擴散層(diffusion layer)和催化層(catalyst layer)，陰陽兩極中間是電解質(zhì)膜(electrolyte membrane)[7-8]。燃料從陽極流道入口流入，氧化劑從陰極流道入口流入，陰陽兩極入口分別在燃料電池的兩個端部，流體相對而行。擴散層和催化層為多孔介質(zhì)區(qū)域，流體擴散到催化層后發(fā)生電化學反應(yīng)，陽極電化學反應(yīng)產(chǎn)生的氫離子通過電解質(zhì)膜進入陰極，電子則通過陽極擴散層及集流板和外電路進入到陰極集流板，再經(jīng)陰極擴散層，進入陰極催化層，在那里和氧化物及氫離子發(fā)生陰極電化學反應(yīng)。陰陽兩極不斷供應(yīng)燃料和氧化劑，電化學反應(yīng)就不斷的進行下去，外電路中的負載就得到了電能的供應(yīng)，這就是簡單的燃料電池將燃料化學能轉(zhuǎn)化為電能的過程。

甲醇燃料電池的陽極、陰極催化層電化學反應(yīng)方程式分別為：

陽極催化層：

CH3OH+H2O=CO2+6e-+6H+

(6)

圖3 PEMFC 結(jié)構(gòu)圖

陰極催化層：

(7)

陽極催化層電化學反應(yīng)產(chǎn)生二氧化碳、氫離子和電子，二氧化碳從陽極出口流出，電子從外電路到達陰極，氫離子通過電解質(zhì)膜到達陰極，這樣就形成一個閉合回路[9]。陰極氧氣和電子及氫離子反應(yīng)產(chǎn)生水，水會反向擴散到陽極。在質(zhì)子交換膜燃料電池中，氫離子的傳輸離不開水，只有膜的納米尺度空隙中存在足夠的水時，水合氫離子才能夠從陽極到達陰極。攜帶氫離子的水為離子拖曳水，這種現(xiàn)象為電滲拖曳[10]。在甲醇燃料電池中，不僅有水的電滲拖曳，還有甲醇的電滲拖曳和甲醇的擴散，到達陰極的甲醇和陰極的氧反應(yīng)，產(chǎn)生寄生電流。

3.2 電化學模型

在PEMFC的陰陽極催化層不斷發(fā)生電化學反應(yīng)，不斷產(chǎn)生和消耗電子和氫離子，反應(yīng)發(fā)生的驅(qū)動力是同一網(wǎng)格上固相電勢(電子電勢)與膜電勢(離子電勢)之間的電勢差，即活化過電勢。電子電勢傳輸主要在固相區(qū)域中，包括集流板、擴散層及催化層，其傳輸方程見式(8)。離子電勢的傳輸發(fā)生在催化層及膜區(qū)域中，其傳輸方程見式(9)。電勢傳輸方程中的擴散系數(shù)即電導率，包括電子在催化層、擴散層、集流板中的電導率，及離子在催化層和膜中的電導率。傳輸方程中的源項R為體電流密度，即單位時間單位體積產(chǎn)生的電量。體電流的大小與反應(yīng)物濃度、活化過電勢、系統(tǒng)溫度、離子交換膜厚度等有關(guān)[11-12]。

·(σsol?sol)+Rsol=0

(8)

(9)

式中，σ為電導率(Ω-1·m-1)，下標sol表示固體(陰陽電極)，m代表膜；?為電勢(V)；R為體電流密度(A·m-3)。

采用簡化的Butler-Volmer 公式(即Tafel 公式)來計算陽極和陰極的交換電流密度分別為：

(10)

(11)

式中，jref為參考交換電流密度(A·m-2)；ζ為反應(yīng)區(qū)比表面積(1/m)；[A]和[A]ref分別為反應(yīng)物濃度及參考濃度 (mol·m-3)；γ為濃度指數(shù)；ηan為過電勢；α為傳遞系數(shù)；F為法拉第常數(shù)(96 485 C/mol)。

在陽極催化層，電子體電流密度等于負的陽極交換電流密度，離子體電流密度等于陽極交換電流密度；在陰極催化層，電子體電流密度等于陰極交換電流密度，離子體電流密度等于負的陰極交換電流密度。陽極陰極的活化過電勢不同，有：

陽極催化層：

Rsol=-Ran(<0)，Rm=+Ran(>0)

陰極催化層：

Rsol=+Rcat(>0) ，Rm=-Rcat(<0)

陽極過電勢：

ηan=?sol-?m

(12)

陰極過電勢：

ηcat=?sol-?m-Vopen

(13)

3.3 質(zhì)量源項

流體混合物組分發(fā)生變化的地方，需要在連續(xù)方程和組分傳輸方程中增加源項。燃料電池陽極催化層和陰極催化層分別發(fā)生電化學反應(yīng)，從反應(yīng)方程式(6)和(7)可以看出，陽極催化層消耗掉甲醇和水，產(chǎn)生二氧化碳，陰極催化層消耗掉氧氣，產(chǎn)生水，當考慮寄生電流時，陰極還會有二氧化碳產(chǎn)生。陽極催化層消耗的甲醇包括電化學反應(yīng)消耗、電滲拖曳和甲醇擴散，其質(zhì)量源項表達式見式(14)；陽極水的變化包括電化學反應(yīng)和電滲拖曳消耗的水及反向擴散來的水，其質(zhì)量源項表達式見式(15)；陽極電化學反應(yīng)產(chǎn)生二氧化碳，其質(zhì)量源項表達式見式(16)。陰極催化層消耗的氧氣包括產(chǎn)生主電流需要的氧氣和產(chǎn)生寄生電流需要的氧氣，其質(zhì)量源項表達式見式(17)；陰極水的變化包括主電流反應(yīng)產(chǎn)生的水、寄生電流反應(yīng)產(chǎn)生的水、電滲拖曳作用帶來的水和反向擴散減少的水，其質(zhì)量源項表達式見式(18)；陰極產(chǎn)生的二氧化碳由于寄生電流而存在，其表達式見式(19)。在連續(xù)方程(mass 方程)中加入變化組分的質(zhì)量源項。

Dch3ohρ(2Y))

(14)

(15)

(16)

(17)

(18)

(19)

式中，S為質(zhì)量源項(kg·m-3·s-1)；M為分子量(kg·m-3) ；ρ為混合物密度(kg·m-3) ；Y為組分質(zhì)量分數(shù)；Rp為寄生電流密度(A·m-3)；ndr為電滲拖曳系數(shù)；Dλ為水含量擴散系數(shù)(mol·m-1·s-1) ；Dch3oh為甲醇擴散系數(shù)(m2·s-1) 。

3.4 能量源項

由于電流的產(chǎn)生，使得燃料電池除流道外的各個區(qū)域的能量方程中都需要增加源項。在催化層，電化學反應(yīng)的焓變并沒有完全轉(zhuǎn)化成電能，為了克服活化能壘有一部分能量損失掉了，催化層能量方程源項見式(20)，在集流板、擴散層和交換膜中，能量方程源項為導電介質(zhì)的散熱量，見式(21)，式中電流密度為單位時間通過單位面積上的電荷，該量為矢量，在某一方向上分量大小為該方向上電勢梯度與電導率的乘積。如果考慮到水蒸氣發(fā)生相變，那么要在擴散層和催化層的能量源項中增加水的汽化潛熱。

Sh=-Δhrect-Ran,catηan,cat

(20)

Sh=I2Rohm

(21)

式中，Δhrect為電化學反應(yīng)焓變(W·m-3)；I為電流密度(A·m-2)，I=σ?；Rohm為電阻，

3.5 自定義內(nèi)存

在不同自定義函數(shù)之間共用的變量值存儲在自定義內(nèi)存中，因此給自定義內(nèi)存賦值位置在求解燃料電池問題中極為中重要，要根據(jù)求解過程中不同宏的調(diào)用順序，在該內(nèi)存值在被引用之前為其賦值。在燃料電池問題中，至少需要如下自定義內(nèi)存：過電勢、交換體電流密度、水含量、膜的電導率、水的電滲拖曳系數(shù)、電流密度及其在x、y、z方向上的分量，以及在計算過程及結(jié)果中感興趣的變量等。過電勢在計算交換體電流密度時使用，在電勢傳輸方程源項中賦值；交換體電流密度在各個方程的源項中都可能用到，所以定義在adjust宏中賦值；水含量在計算膜電導率、水含量擴散系數(shù)及水的電滲拖曳系數(shù)時都會用到，所以在混合物屬性宏中賦值；膜的電導率及水的電滲拖曳系數(shù)也在屬性宏中賦值；電流密度及其分量只在能量方程源項中引用，且除流道域外各個區(qū)域的能量方程均勻源項，所以在能量方程源項中給電流密度賦值即可。在計算過程陽極陰極的電流不斷接近，收斂時達到相等，我們既關(guān)心電流過程值便于調(diào)試程序，也關(guān)心計算結(jié)果電流的大小，所以在自定義內(nèi)存中要有電流，電流是催化層所有網(wǎng)格體電流密度的總和，在adjust宏中給電流賦值。

4 程序調(diào)試

電化學反應(yīng)主要通過電勢標量傳輸方程及其他方程增加源項體現(xiàn)出來的。程序調(diào)試時，首先不考慮電化學反應(yīng)，即不求解電勢標量傳輸方程也不增加源項的情況下，確?；旌衔锪鲌龅恼_性。要保證陽極入口流入的流體在陽極流道、陰極入口流入的流體在陰極流道，一是要設(shè)定膜域中由于壓力差引起的宏觀流動的速度為零，二是設(shè)定合理的混合物擴散系數(shù)。待無電化學反應(yīng)的流場計算正確后，再加載adjust 宏及電勢傳輸方程中源項及質(zhì)量源項，調(diào)試電流密度場，最后加載能量方程中源項，調(diào)試溫度場。

5 計算結(jié)果

開路電壓為0.95 V，進料溫度353 K，陰陽極參考電流密度為1 000 A/m2，操作壓力為101 325 Pa，陽極入口流體質(zhì)量流量為5×10-6kg/s，陽極入口流體質(zhì)量流量為6×10-7kg/s，陰陽極電荷傳遞系數(shù)為2，電池擴散層、催化層及膜的氣孔率為0.5，陽極濃度指數(shù)為0.5，陰極濃度指數(shù)為1。流道為單通道時計算得到的電流密度場見圖4、溫度場見圖5、甲醇質(zhì)量場見圖6、水分質(zhì)量場見圖7。

圖4 電流密度分布

圖5 溫度分布

圖6 甲醇質(zhì)量組分分布

圖7 水分質(zhì)量組分分布

6 總結(jié)

用Fluent 流體計算軟件求解質(zhì)子交換膜燃料電池問題，關(guān)鍵在于設(shè)定正確的自定義電勢標量擴散系數(shù)(即電導率)、標量傳輸方程源項(體電流密度)、連續(xù)方程及組分傳輸方程中質(zhì)量源項、能量方程中能量源項，還要設(shè)定合理的自定義內(nèi)存變量及合適的賦值位置。

[1] 奧黑爾.燃料電池基礎(chǔ)[M].王曉紅，黃宏，譯.北京：電子工業(yè)出版社，2007.

[2] 韓占鐘，王敬，蘭小平.Fluent流體工程仿真計算實例與應(yīng)用[M].北京：北京理工大學出版社，2004.

[3] 王瑞金，張凱，王剛.Fluent技術(shù)基礎(chǔ)與應(yīng)用實例[M].北京：清華大學出版社，2007.

[4] 張建文，楊振亞，張政.流體流動與傳熱過程的數(shù)值模擬基礎(chǔ)與應(yīng)用[M].北京：化學工業(yè)出版社，2009.

[5] Sukkee U m,Wang C Y,Chen K S.Computational Fluid Dynamics Modelingof Proton Exchange Membrane Fuel Cells[J].J Electrochemical Society,2000，147(12):4485-4493.

[6] Springer T E,Zawodzinski T A,Gottesfeld S.Polymer Electolyte Fuel CellModel[J].J Electrochemical Soc,1991，138(8):2334-2342.

[7] 熊濟時，肖金生，潘牧.質(zhì)子交換膜燃料電池流場模擬與結(jié)構(gòu)尺寸優(yōu)化[J].武漢理工大學學報：交通科學與工程版,2009,33(6)：534-537.

[8] 馬利軍，林才順，薛方勤，等.質(zhì)子交換膜燃料電池流場截面設(shè)計及三維模擬[J].濕法冶金，2008,27(1):52-55.

[9] 劉桂成，張浩，王一拓，等.液體進料直接甲醇燃料電池工況的匹配及優(yōu)化[J].電池，2012,42(1):7-10.

[10] 潘科，齊亮，劉洋，等.質(zhì)子交換膜的電滲系數(shù)[J].電池，2008,38(2):75-78.

[11] Kulikovsky A A,Divisek J,Kornyshev A A.Modeling the Cathode Com-partment of Polymer Electrolyte Fuel Cells: Dead and Active Reaction Zones[J].J Electrochemical Society,1999，146(11):3981-3991.

[12] Mazumder S，Cole J V.Rigorous 3-D Mathematical Modeling of PEM Fuel Cells II.Model Predictions with Liquid Water Transport[J].J Electrochemical Soc,2003，150(11):A1510-A1517.