翟占虎，李 巍，劉國(guó)棟

(1.哈爾濱工業(yè)大學(xué)能源科學(xué)與工程學(xué)院，黑龍江 哈爾濱 150001；2.海洋石油工程股份有限公司，天津 300452)

太陽(yáng)能因其儲(chǔ)量豐富、清潔可再生等特性，成為當(dāng)前研究的熱點(diǎn)問(wèn)題，但是由于其熱流密度低、具有周期性且不易儲(chǔ)存，目前還只是作為輔助能源.近年來(lái)太陽(yáng)能聚焦技術(shù)取得長(zhǎng)足發(fā)展，將光熱轉(zhuǎn)換和高溫?zé)峄瘜W(xué)反應(yīng)相耦合的光化學(xué)利用技術(shù)引起重視，該技術(shù)以氧載體材料為依托，在高溫下將 CO2、H2O 等熱分解產(chǎn)生H2和CO等碳?xì)淙剂?

圍繞著太陽(yáng)能反應(yīng)器，國(guó)內(nèi)外學(xué)者開(kāi)展大量的研究工作.瑞士PSI太陽(yáng)能技術(shù)實(shí)驗(yàn)室Steinfeld[1]教授領(lǐng)導(dǎo)的團(tuán)隊(duì)設(shè)計(jì)并試驗(yàn)的太陽(yáng)能熱化學(xué)反應(yīng)器，能夠?qū)H4-ZnO在高溫下分解為Zn和合成氣.Chueh等[2]設(shè)計(jì)的太陽(yáng)能反應(yīng)器，參考植物將太陽(yáng)能轉(zhuǎn)化為燃料的過(guò)程，研究了CeO2在熱化學(xué)循環(huán)中生產(chǎn)燃料的熱力學(xué)和動(dòng)力學(xué)性質(zhì)以及反應(yīng)的穩(wěn)定性，在不進(jìn)行廢熱回收的前提下，該反應(yīng)的熱化學(xué)學(xué)循環(huán)效率理論上可達(dá)到16%～19%.Meier等[3]設(shè)計(jì)并實(shí)驗(yàn)的太陽(yáng)能間接加熱回轉(zhuǎn)窯，可以在1 200 K～1 400 K溫度區(qū)間內(nèi)實(shí)現(xiàn)高溫煅燒反應(yīng)，并獲得高純度石灰.吳璽[4]針對(duì)ZnO與CH4-Ar的重整制氫反應(yīng)，數(shù)值模擬研究了吸收腔內(nèi)溫度場(chǎng)分布情況，并分析了腔體溫度、反應(yīng)物料進(jìn)口速度和CH4-Ar體積比對(duì)重整制氫反應(yīng)效率的影響.王富強(qiáng)等[5]采用數(shù)值模擬對(duì)比研究了輻射模型多孔介質(zhì)太陽(yáng)能接收器傳熱性能的影響.黃興等[6]研究了催化劑NiFe2O4的熱解機(jī)理并應(yīng)用于太陽(yáng)能熱化學(xué)反應(yīng)器，采用Ansys Fluent研究了反應(yīng)器內(nèi)速度場(chǎng)和溫度場(chǎng)，分析了顆粒直徑和氣體流速對(duì)反應(yīng)器內(nèi)部的溫度分布影響.Lougou等[7]數(shù)值模擬研究了太陽(yáng)能熱化學(xué)反應(yīng)器作用下的甲烷重整反應(yīng)，分析了NiFe2O4催化作用下，混合氣體流速、操作壓力、催化劑裝載量和輻射溫度分布等參數(shù)對(duì)反應(yīng)轉(zhuǎn)化率的影響.張浩等[8]采用數(shù)值模擬結(jié)合反應(yīng)動(dòng)力學(xué)，研究了太陽(yáng)能反應(yīng)器中的甲烷重整反應(yīng)和積碳特性，得到多孔介質(zhì)的最佳空隙率和顆粒直徑.Bader等[9]設(shè)計(jì)了太陽(yáng)能高溫?zé)峄瘜W(xué)反應(yīng)器，以二氧化鈰(CeO2)為固定床催化劑，利用該金屬氧化還原作用的分解CO2，該反應(yīng)器中催化材料與吸收腔不直接接觸，CeO2相對(duì)于Fe3O4和ZnO等催化材料，具有儲(chǔ)放氧性能高、比表面積大、傳熱傳質(zhì)速率高和結(jié)構(gòu)穩(wěn)定等優(yōu)點(diǎn)[10]，避免了催化劑的污染和分離回收，同時(shí)獨(dú)立的反應(yīng)器提升了系統(tǒng)穩(wěn)定性和操作性，適用于存在一定溫度波動(dòng)的高溫反應(yīng)器.

太陽(yáng)能熱化學(xué)反應(yīng)器反應(yīng)效率受熱流密度分布和金屬氧化物催化劑的氧化還原過(guò)程的影響，實(shí)驗(yàn)直接獲得相關(guān)的參數(shù)比較困難，因此本文數(shù)值模擬研究Bader所設(shè)計(jì)的太陽(yáng)能熱化學(xué)反應(yīng)器，反應(yīng)過(guò)程中氧空位濃度δ的變化通過(guò)求解編寫的UDS(用戶自定義標(biāo)量方程)實(shí)現(xiàn)，質(zhì)量和能量變化通過(guò)UDF(用戶自定義函數(shù))加載到多孔介質(zhì)區(qū)域，最后分析氧化還原反應(yīng)階段的溫度、氧空位以及產(chǎn)物濃度分布情況.

1 反應(yīng)器結(jié)構(gòu)與數(shù)學(xué)模型

1.1 反應(yīng)器結(jié)構(gòu)及物理模型

反應(yīng)器整體結(jié)構(gòu)如圖1所示，由一石英窗、光圈、保溫層、圓筒形吸收腔及腔內(nèi)環(huán)形排列著的七個(gè)反應(yīng)器構(gòu)成，吸收腔內(nèi)部長(zhǎng)為347 mm，半徑為152.5 mm，每個(gè)反應(yīng)器由同軸圓管構(gòu)成，反應(yīng)器外管的外/內(nèi)半徑分別為34.9 mm/31.8 mm，內(nèi)管的外/內(nèi)半徑分別為22.2 mm/19.1 mm，圓筒壁面厚為3.1 mm，兩圓筒夾層厚為9.5 mm，夾層中放置718 g固定的金屬氧化物催化劑二氧化鈰(CeO2)，該催化劑為顆粒直徑為5 mm，催化劑顆粒的孔隙率為70%，固定床層的空隙率為40%.

圖1 反應(yīng)器幾何結(jié)構(gòu)圖

反應(yīng)器工作流程如下，經(jīng)過(guò)聚焦后的高通量太陽(yáng)輻射能(3 kW)，透過(guò)光圈射在吸收腔壁面和七個(gè)圓管反應(yīng)器的表面，吸收了太陽(yáng)輻射能的氧化鋁管壁將熱量傳入固定床層的二氧化鈰，催化材料二氧化鈰(CeO2)其表面及體相的晶格氧原子能夠直接參與反應(yīng)并被消耗[11]，同時(shí)形成氧空位(δ)，產(chǎn)出化學(xué)燃料.在反應(yīng)器工作過(guò)程中，周期性循環(huán)地向反應(yīng)管內(nèi)通入N2和CO2，通入N2時(shí)長(zhǎng)為100 s，質(zhì)量流量為1.1 g·s-1，該階段發(fā)生還原反應(yīng)，出口產(chǎn)物有N2和O2；通入CO2時(shí)長(zhǎng)為100 s，質(zhì)量流量為1.7 g·s-1，該階段發(fā)生氧化反應(yīng)，出口產(chǎn)物有CO2和CO.反應(yīng)方程式如下所示：

還原反應(yīng)(Reduction)

(1)

氧化反應(yīng)(Oxidation)

CeO2-δrd+(δrd-δox)CO2=CeO2-δox+(δrd-δox)CO

(2)

二氧化鈰顆粒及其所構(gòu)成的固定床，在模擬中被簡(jiǎn)化為多孔介質(zhì)，根據(jù)相關(guān)公式計(jì)算出的輸運(yùn)性質(zhì)如表1所示.

表1 二氧化鈰和固定床層輸運(yùn)性質(zhì)

為簡(jiǎn)化計(jì)算，取模型的七分之一進(jìn)行模擬，首先模擬得到吸收腔內(nèi)熱流密度分布，然后將熱流密度加載到單個(gè)反應(yīng)器，再進(jìn)行反應(yīng)器內(nèi)流動(dòng)、傳熱和化學(xué)反應(yīng)模擬.

1.2 控制方程

反應(yīng)過(guò)程中滿足連續(xù)性方程、動(dòng)量方程、組分輸運(yùn)方程、能量守恒方程.

連續(xù)性方程

(1)

公式中：ρ、u、Ri分別為密度、速度和質(zhì)量源項(xiàng)；下標(biāo)f、s、i分別為流體、固體和組分i.

動(dòng)量守恒方程

(2)

公式中：φ、p、μ、K、CE分別為空隙率、壓力、動(dòng)力粘度、滲透率和慣性阻力系數(shù).

組分輸運(yùn)方程

(3)

公式中：φeff、Yi分別為有效空隙率和組分i的質(zhì)量分?jǐn)?shù).

能量守恒方程

(4)

(5)

(6)

DO輻射方程：

(7)

由于反應(yīng)器吸收腔內(nèi)溫度較高，需要考慮其內(nèi)部的輻射特性，在模擬吸收腔內(nèi)部熱流密度分布時(shí)采用DO輻射模型.

1.3 化學(xué)反應(yīng)動(dòng)力學(xué)方程

目前的文獻(xiàn)中，沒(méi)有準(zhǔn)確描述二氧化鈰催化分解CO2的反應(yīng)動(dòng)力學(xué)表達(dá)式，本研究中采用的動(dòng)力學(xué)表達(dá)式主要參考Bulfie和Rohini Bala的文章[14-15].在二氧化鈰參與氧化還原反應(yīng)的過(guò)程中，其非化學(xué)計(jì)量數(shù)δ決定反應(yīng)進(jìn)行速率的大小，以下為δ隨時(shí)間變化的微分方程，以及各組分質(zhì)量源表達(dá)式.

δ與t的關(guān)系：

(8)

(9)

在氧化與還原過(guò)程中，CeO2既能釋放O2，也能吸收O2，因此氧氣質(zhì)量源項(xiàng)的表達(dá)式如下：

(10)

研究表明[17]，高溫下CO2的直接分解作用也很明顯，且為可逆反應(yīng)，各組分質(zhì)量源項(xiàng)的表達(dá)式如下：

(11)

公式中：i=O2，CO，CO2，相應(yīng)地vi=0.5，1，-1.

O2的凈產(chǎn)量為CO2分解所釋放的量與CeO2氧化消耗的量之和，如下所示：

RO2，net=RO2，diss+RO2，ceria

(12)

預(yù)計(jì)熱分解產(chǎn)物CO會(huì)在短時(shí)間內(nèi)和O2重新結(jié)合成CO2，故本研究認(rèn)為氧化階段產(chǎn)生的O2，一部分被CeO2吸收來(lái)填補(bǔ)氧空位，其余的全部與CO重新結(jié)合.

本研究以上述的非化學(xué)計(jì)量數(shù)δ的微分方程、各組分的質(zhì)量源項(xiàng)表達(dá)式為理論基礎(chǔ)，來(lái)實(shí)現(xiàn)上述化學(xué)反應(yīng).

2 模擬結(jié)果分析

Bader[9]文中給出的太陽(yáng)能熱化學(xué)反應(yīng)器為概念設(shè)計(jì)，并未進(jìn)行實(shí)驗(yàn)測(cè)量，反應(yīng)器在設(shè)計(jì)運(yùn)行參數(shù)下，工作溫度約為1 800 K，CO產(chǎn)率在預(yù)估為3.9×10-4mol·s-1，模擬得到的溫度分布如圖4和圖11所示，以及出口CO摩爾產(chǎn)率如圖17所示，較為接近，可認(rèn)為模擬所采用模型是合理的.

2.1 吸收腔內(nèi)熱流密度分布

吸收腔內(nèi)部熱流密度分布情況，如圖2所示.可以看出圓管反應(yīng)器壁面的熱流量明顯高于吸收腔壁面，經(jīng)計(jì)算圓管反應(yīng)器壁面的平均熱流密度為4 600 W·m-2，總熱流量為318.5 W，吸收腔壁面的熱流密度為181.5 W，熱量損失占36.3%.沿著圓管反應(yīng)器軸的熱流密度分布曲線，如圖3所示.沿著軸向先增高后減少，并且在x=0.12 m達(dá)到最大值6 900 W·m-2.

圖2 反應(yīng)器熱流密度分布云圖圖3 沿X軸向熱流密度

2.2 還原階段分析

(1)反應(yīng)器溫度分析

t=1 s和100 s時(shí)反應(yīng)器溫度分布云圖，如圖4所示.從左側(cè)圖可看出反應(yīng)器壁面溫度最高為1 880 K，t=100 s時(shí)超過(guò)1 860 K溫度的范圍明顯少于t=1 s時(shí)的，這是因?yàn)镃eO2的還原反應(yīng)是吸熱反應(yīng)，隨著反應(yīng)的經(jīng)行，反應(yīng)器溫度下降.

圖4 t=1 s和100 s時(shí)反應(yīng)器溫度分布

反應(yīng)器壁面軸向溫度在不同時(shí)刻的溫度分布和多孔介質(zhì)的體積平均溫度隨時(shí)間的變化，如圖5、圖6所示.壁面溫度和體積平均溫度都隨著反應(yīng)的進(jìn)行，不斷減小，進(jìn)一步驗(yàn)證了CeO2的還原反應(yīng)是吸熱反應(yīng).

圖5 反應(yīng)器壁面軸向溫度在不同時(shí)刻的溫度分布圖6 多孔介質(zhì)體積平均溫度隨時(shí)間的變化

(2)氧空位濃度δ與O2產(chǎn)率

t=1 s時(shí)和t=100 s多孔介質(zhì)區(qū)域的δ分布云圖，如圖7所示.在t=100 s時(shí)δ最高到0.041，且最高位置在出現(xiàn)在剛進(jìn)入多孔介質(zhì)區(qū)域的地方.這是因?yàn)闅怏w剛進(jìn)入多孔介質(zhì)時(shí)，由于不含氧氣，故此處氧氣的分壓為零，遠(yuǎn)低于平衡壓力，CeO2中的氧原子在高溫條件下能夠以氧氣的形式釋放出來(lái)，故此處的氧空位最高.從圖7中還可看出，靠近出口處多孔介質(zhì)的δ最小，這是因?yàn)檠鯕獾牟粩噌尫牛鯕庠谙掠蔚馁|(zhì)量分?jǐn)?shù)達(dá)到最高，即氧氣分壓也達(dá)到最高，故此處的CeO2中的氧原子受氧氣分壓的抑制，釋放速率較慢.

圖7 不同時(shí)刻δ在多孔介質(zhì)區(qū)域的分布云圖

沿軸向δ是不斷減小的，且隨著時(shí)間增加，δ的增長(zhǎng)變緩，多孔介質(zhì)上游區(qū)域的δ增長(zhǎng)率明顯高于下游，如圖8所示.由圖9可知，氧氣的質(zhì)量分?jǐn)?shù)沿軸向不斷變大，隨之時(shí)間增加，氧氣的質(zhì)量分?jǐn)?shù)不斷降低.這是因?yàn)楫?dāng)CeO2剛開(kāi)始反應(yīng)時(shí)，由于多孔介質(zhì)區(qū)域各處的氧氣分壓均為零，O2的產(chǎn)率很高，隨著反應(yīng)的進(jìn)行，氧氣分壓升高，特別是多孔介質(zhì)下游區(qū)域，其氧氣分壓明顯高于上游，更接近平衡壓力，當(dāng)氧氣分壓越接近平衡壓力，δ增長(zhǎng)越慢，氧氣的產(chǎn)率越小，相反地上游區(qū)域氧氣分壓始終較低，δ增長(zhǎng)始終較快，所以x=0 m處的δ值明顯高于x=0.374 m處的.

圖8 軸向氧空位濃度δ分布(不同時(shí)刻)圖9 軸向氧氣質(zhì)量分?jǐn)?shù)分布(不同時(shí)刻)

出口氧氣的摩爾流率，如圖10所示.在剛通入N2的時(shí)候，出口氧氣的摩爾流率可達(dá)到9.42×10-5mol·s-1，且隨著反應(yīng)時(shí)間的增加，氧氣流率逐漸下降，這是因?yàn)榭傮w上氧氣的分壓越來(lái)越接近平衡壓力(大約為240 Pa)，δ增長(zhǎng)減緩，從CeO2中釋放出的氧氣變少.

圖10 出口氧氣的摩爾流率

2.3 氧化階段分析

(1)反應(yīng)器內(nèi)溫度分布

t=101 s和200 s時(shí)反應(yīng)器溫度分布云圖，如圖11所示.在整個(gè)時(shí)間段內(nèi)，反應(yīng)器的最高溫度為1 770 K，低于還原階段的1 880 K，而從101 s到200 s這段時(shí)間，反應(yīng)器內(nèi)高溫區(qū)域變大.這是因?yàn)橥ㄈ隒O2后，CO2質(zhì)量流量比N2的質(zhì)量流量大，且CO2氣體的比熱容大于N2的比熱容，故在一定的熱流密度下，氧化階段的流場(chǎng)溫度低于還原階段時(shí)的溫度.此外，在整個(gè)氧化階段，反應(yīng)是放熱的，所以隨著時(shí)間的增加，反應(yīng)器整體升溫.

圖11 不同時(shí)刻壁面和計(jì)算域內(nèi)溫度分布云圖

反應(yīng)器壁面軸向溫度在不同時(shí)刻的溫度分布，如圖12所示.從圖12可知，隨著軸向距離的增加，壁面溫度先增高后下降.因?yàn)闅怏w被壁面不斷加熱，隨著流動(dòng)距離溫度變高，在臨近出口時(shí)，由于與入口低溫氣體(1 625 K)換熱，所以溫度出現(xiàn)了下降，由于反應(yīng)放熱，所以隨著時(shí)間增加，壁面溫度升高.從圖13可知，隨著時(shí)間的增加，多孔介質(zhì)體積平均溫度不斷增加，進(jìn)一步說(shuō)明了氧化階段發(fā)生的反應(yīng)是放熱的.

圖12 反應(yīng)器壁面軸向溫度分布圖13 多孔介質(zhì)體積平均溫度隨時(shí)間的變化

(2)氧空位濃度δ與CO產(chǎn)率

從圖14可看出，隨著軸向距離增加，氧空位濃度δ逐漸降低，且隨著時(shí)間的增加，δ也在減少，特別是x=0 m處，δ下降最快，在t=200 s時(shí)多孔介質(zhì)區(qū)域δ的分布云圖，如圖15所示.可見(jiàn)在x=0 m處，δ約為0.03，體積平均δ約為0.020，與還原階段的初始值較接近，這是因?yàn)檠趸A段，CO2高溫下發(fā)生熱分解，釋放出大量O2，使得此時(shí)的氧氣分壓大于平衡壓力，CeO2將O2吸收來(lái)填補(bǔ)自身的氧空位，所以δ出現(xiàn)下降趨勢(shì)，且δ越大，氧氣平衡壓力就越小，CeO2對(duì)O2的吸收速率更高，因此在x=0 m處δ下降最快.

圖14 軸向氧空位濃度δ分布情況(不同時(shí)刻)圖15 t=200s時(shí)δ分布云圖

從圖16可看出，隨著x的增加，CO質(zhì)量分?jǐn)?shù)先增加后減少，且隨著時(shí)間的增加，CO質(zhì)量分?jǐn)?shù)不斷下降，這是因?yàn)楫?dāng)通入CO2時(shí)，CO2熱分解釋放的O2被二氧化鈰CeO2吸收用以填補(bǔ)氧空位，進(jìn)而促進(jìn)CO2的熱分解，所以從x=0 m到x=0.15 m處，CO質(zhì)量分?jǐn)?shù)不斷升高；而CO質(zhì)量分?jǐn)?shù)在x=0.15 m之后出現(xiàn)下降趨勢(shì)，是因?yàn)檫@段區(qū)域的氧空位濃度δ較小，對(duì)O2的吸收能力有限，導(dǎo)致CO2熱分解產(chǎn)生的O2沒(méi)有被完全吸收，與CO再結(jié)合，導(dǎo)致CO質(zhì)量分?jǐn)?shù)下降.

CO摩爾流率隨時(shí)間的變化曲線，如圖17所示.隨時(shí)間增加，CO摩爾流率逐漸下降，這是因?yàn)殡S著反應(yīng)的進(jìn)行，氧空位濃度δ逐漸減少，對(duì)CO2熱分解反應(yīng)的促進(jìn)表少，導(dǎo)致CO產(chǎn)率逐漸減少.

圖16 軸向CO質(zhì)量分?jǐn)?shù)分布 圖17 出口CO的摩爾流率

3 結(jié) 論

本文以Bader設(shè)計(jì)的反應(yīng)器為基礎(chǔ)，建立了用于二氧化鈰小溫差下氧化還原循環(huán)分解CO2的太陽(yáng)能高溫反應(yīng)器的三維瞬態(tài)模型，考慮了流動(dòng)過(guò)程和非化學(xué)計(jì)量反應(yīng)的相互作用及其對(duì)溫度、壓力、組分濃度和燃料產(chǎn)率的影響.研究發(fā)現(xiàn)，反應(yīng)器熱流密度集中分布在靠近光圈的一側(cè)，熱量損失取決于吸收腔壁面溫度和材料保溫性能，應(yīng)盡量減少輻射能量直接作用于吸收腔壁面；反應(yīng)器溫度與反應(yīng)階段密切相關(guān)，還原階段溫度上升，氧化階段溫度下降；還原階段，二氧化鈰釋放氧氣，氧空位濃度δ受氧氣的抑制作用，在氣體入口處最大，出口處最??；氧化階段，二氧化碳在高溫下熱分解，二氧化鈰吸收分解后的氧原子，促進(jìn)了CO2的熱分解作用；氧氣釋放量和二氧化碳熱分解作用與氧空位濃度、氣體流率、反應(yīng)時(shí)間和溫度密切相關(guān)，在太陽(yáng)能反應(yīng)器設(shè)計(jì)階段需綜合考慮上述因素.