張國艷,安興才,陳作雁,劉 剛
(甘肅自然能源研究所,甘肅 蘭州 730046)
光催化技術在治理含有機染料的廢水行業(yè)具有廣泛的應用前景[1-3],將光催化劑固定在玻璃珠上,既能保證較好的光催化活性、透光性,又解決了催化劑的回收問題,而且能夠多次使用、連續(xù)處理廢水[4-5]。負載型TiO2玻璃珠組成的固定床反應器在光催化領域受到極大關注,其中固定床內流場和壓力場是評價固定床設備的重要參數(shù),二者直接影響反應器的效率,所以能夠準確地預測固定床內的流場和壓力場的變化規(guī)律對固定床的設計和操作有重要的意義。
目前國內外許多學者采用計算流體動力學(CFD,computational fluid dynamics)的方法對隨機或規(guī)則填充固定床的流場和壓力場進行了模擬計算,重點考察了管徑比、雷諾數(shù)、孔隙率、顆粒形狀、壁面效應等對流場和壓力場的影響[6-14]。這些研究均以部分填充床為基礎,不能夠真實全面地反應整根管的流場和壓力場的情況,并且固定床形狀為圓柱形或直方形,對環(huán)形固定床的流場和壓力場研究仍然較少。環(huán)形固定床相對于圓柱形或直方形固定床,內部光源分布均勻、光效率高。然而,環(huán)形固定床反應器的研究基本以實驗為主,通過實驗手段雖然能直接獲得床層內壓降,但其內部的壓降和流場的具體變化規(guī)律無法形象和實時地表達,且環(huán)形固定床內質量傳遞不明確。通過數(shù)值模擬的方法模擬整根反應管的流場和壓力場,不僅可以觀測到床層內部流場和壓降變化、質量傳遞情況,還可以優(yōu)化實驗工況參數(shù)。
首先研究環(huán)形固定床反應器的內部阻抗值與玻璃珠直徑、入口流速的關系,然后運用CFD方法對環(huán)形固定床反應器內部的流場和壓降進行數(shù)值計算,研究玻璃珠直徑和入口流速對流場和壓力場的影響,并對模擬計算結果進行了圖譜分析,結合實驗確定出最佳玻璃珠直徑和入口流速,優(yōu)化了反應器結構,為反應動力學模擬及工程應用提供技術基礎。
環(huán)形固定床光催化氧化反應器基本單元和實驗裝置見圖1。該反應器由外管、內管、內外管間填充載有光催化劑的玻璃珠、下端面入口和上端面出口等幾部分構成,其結構參數(shù)如表1所列。在實驗裝置中反應器與水平面呈36°夾角放置,實驗光源位于內外管中心線處,其長度為1.5 m。
圖1 環(huán)形固定床光催化反應器基本單元示意圖和實驗裝置圖Fig.1 Schematic diagram of basic unit and experimental device layout of annular fixed bed photocatalytic reactor
表1 環(huán)形固定床各部分名稱及尺寸
由于光催化過程中光傳遞效率的限制,環(huán)形固定床內玻璃珠層數(shù)不宜過多,在以上結構參數(shù)條件下,利用CFD軟件模擬研究徑向1~2層玻璃珠的環(huán)形固定床反應器的流速場和壓力場,每次裝填相同直徑玻璃珠,其玻璃珠直徑分別為4.88 mm、5.86 mm、6.93 mm、7.91 mm。
為了更好地描述環(huán)形固定床內顆粒的分布對流體流場和壓力場的影響,考慮到固定床的形狀,采用軸對稱數(shù)學模型。連續(xù)性方程、動量方程、湍動能方程和湍動能耗散效率方程在軸對稱穩(wěn)態(tài)情況下,通用形式[15]為
(1)
其中:Γφ為φ的擴散系數(shù),對不同的方程表示不同的值;Sφ為源項,包括壓力梯度項,在固定床內表示流體與固定顆粒摩擦的阻力項。
對圓柱形或直方形固定床的流場和壓力場的研究[9-13,16-17],主要考慮固定床顆粒直徑、顆粒比表面積、顆粒形狀等的影響。Ergun[18]提出計算床層內壓降的公式為
(2)
其中:ΔP為床層壓降(Pa);L為床層高度(m);ρ為流體密度(kg/m3);Us為入口流速(m/s);μ為動力粘度(Pa·s);dp為顆粒直徑(m);ε為床層孔隙率,無量綱。研究利用式(2)計算環(huán)形固定床反應器的壓降。
環(huán)形管內部填充玻璃珠,玻璃珠之間、玻璃珠與壁面之間彼此相連,為了保證網(wǎng)格質量,將玻璃珠直徑縮小1%。這樣處理會導致玻璃珠填充率的誤差為3%,壓降誤差為10%,均在工程應用誤差范圍之內[13]。
應用ICEM軟件進行建模并劃分網(wǎng)格,采用非結構四面體網(wǎng)格。環(huán)形固定床端面和L/2斷面網(wǎng)格示意圖如圖2所示,由圖2可以看到端面處網(wǎng)格較均勻,而L/2斷面的網(wǎng)格扭曲率(skewness)≤0.7,網(wǎng)格質量(quality)≥0.3,這是因為玻璃珠之間、玻璃珠與壁面之間間隙僅為0.03 mm,出現(xiàn)網(wǎng)格扭曲率較大、網(wǎng)格質量較差的情況是不可避免的。由文獻[10,19]可知,skewness≤0.7和quality≥0.3可以保證網(wǎng)格質量,且數(shù)值計算結果可靠。
(1) 網(wǎng)格無關性驗證 以直徑為5.86 mm的玻璃珠為例,驗證網(wǎng)格無關性和時間步長獨立性。
圖2 環(huán)形固定床端面和L/2斷面網(wǎng)格示意圖Fig.2 End face of the annular fixed bed and L/2section grid diagram
研究模擬關注的對象為流速,選取6個點比較具體的流速,其坐標分別為1(22,10,0)、2(22,10,300)、3(22,10,600)、4(22,10,900)、5(22,10,1250)、6(22,10,1500)。4種計算網(wǎng)格數(shù)見表2。對表2中的4種網(wǎng)格數(shù)進行數(shù)值模擬,6個點的模擬結果見圖3。
表2 4種計算網(wǎng)格數(shù)
從圖3可以看出,隨著網(wǎng)格數(shù)的增加,6個點的速度在網(wǎng)格數(shù)為1 441萬時的波動比較大,其他3種情況對計算結果影響很小,為保證計算結果的準確性和計算效率,取2 536萬的網(wǎng)格作為計算網(wǎng)格。
(2) 時間步長獨立性驗證 時間步長分別取0.025 s、0.05 s、0.1 s、0.2 s進行數(shù)值模擬,6個點的速度隨時間步長的變化如圖4所示。
從圖4可以看出,6個點在不同的時間步長下結果基本不變,綜合考慮計算周期和計算精度,計算的時間步長取0.1 s。
圖3 網(wǎng)格無關性驗證Fig.3 Grid independence verification
圖4 時間步長獨立性驗證Fig.4 Time step independence verification
研究中進行的是不可壓縮液體單相流的模擬計算,水從反應器下端面入口進入,經(jīng)過玻璃珠固定床,到達上端面出口。入口條件為速度入口;出口條件為壓力出口,大小為0,即為大氣壓;玻璃珠及邊壁均采用無滑移條件,各項流速為0。
采用基于壓力的求解器,SIMPLE算法,二階迎風差分格式,收斂標準為5×10-4。
由于環(huán)形固定床是由玻璃珠堆積而形成的多孔結構,在數(shù)值模擬過程中需要開啟Porous Zone選項,對多孔介質區(qū)域進行設置。在Relative Velocity Resistance Formulation中設置為0;在Inertial Resistance中X、Y、Z方向分別輸入內部阻抗值,不同直徑玻璃珠在不同入口流速時的內部阻抗值如圖5所示,其他保持默認設置。
圖5 不同直徑玻璃珠的內部阻抗值隨入口流速的變化Fig.5 Variation of internal impedance of glass beads withdifferent diameters with inlet flow rate
(1) 內部阻抗值的確定 應用Ergun方程對環(huán)形固定床的壓降進行了理論計算,根據(jù)理論計算值確定內部阻抗值與不同玻璃珠直徑和入口流速的關系,其結果如圖5所示。
由圖5可以看出,不同直徑玻璃珠的內部阻抗值隨入口流速的增大而逐漸增加,當流速小于0.20 m/s時,內部阻抗值隨入口速度的增加呈拋物線趨勢增加;流速大于0.20 m/s時,內部阻抗值趨于穩(wěn)定;同一入口流速下,隨著玻璃珠直徑的增加,內部阻抗值呈階梯狀減小。不同玻璃珠直徑的內部阻抗值與入口速度的擬合關系為
y=y0+(A/(w×sqrt(pi/2)))×
exp(-2×((v-vc)/w)2)。
(3)
式(3)中各參數(shù)取值見表3。
表3 擬合關系式的參數(shù)
(2) 內部阻抗值的驗證 對內部阻抗值與不同直徑玻璃珠和入口流速的關系進行了數(shù)值模擬,并做了相應的實驗,其結果如圖6所示。不同直徑玻璃珠的數(shù)值模擬、實驗和理論壓差隨流速的變化趨勢基本相同,均呈拋物線趨勢增加。入口流速相同時,玻璃珠直徑越大,壓差越小。流速大于0.15 m/s時,模擬壓差、實驗壓差與理論壓差基本相等;當流速小于0.15 m/s時,三者最大誤差為10%。因此,內部阻抗值可以保證數(shù)值模擬結果的準確性。
圖6 不同直徑玻璃珠的實驗壓差、理論壓差、模擬壓差隨流速的變化Fig.6 Experimental/theoretical/simulated pressure differenceof glass beads with different diameters with flow rate
雷諾數(shù)判斷環(huán)形固定床內流體的流動狀態(tài):
其中:v為流體斷面平均流速(m/s);d為固定床當量直徑(m);υ為運動黏度(m2/s)。經(jīng)計算,當v=0.115 m/s時,Re=2 300,管內流體達到紊流狀態(tài)。
研究以直徑為5.86 mm的玻璃珠為例,考察不同入口流速對流場和壓力場的影響。
(1) 入口流速對流場的影響 不同入口流速L/2斷面處流速云圖分布如圖7所示,隨著入口流速的增加,顆粒間隙內流速逐漸增加,并且隨著入口流速的增加,環(huán)形管內流體流速和玻璃珠表面速度梯度逐漸增大。入口速度分別為0.10 m/s、0.15 m/s、0.20 m/s、0.22 m/s時,管內主流速大小分別為0.20 m/s、0.35 m/s、0.45 m/s、0.45 m/s,且入口速度由0.15 m/s變化到0.20 m/s時,管內流速變化顯著,玻璃珠表面速度梯度明顯增大,因此,可認為入口流速為0.15~0.20 m/s是玻璃珠表面速度梯度增大的轉折區(qū)間。
環(huán)形固定床環(huán)間距10 mm,管長1.50 m,內部填充玻璃珠直徑為5.86 mm時,不同入口流速的環(huán)形管內流體流速分布如圖8所示,隨著入口流速的增加,管內流體流速分布比例向流速較大的方向移動且趨于緩和。入口流速分別為0.10 m/s、0.15 m/s、0.20 m/s、0.22 m/s時,流體流速≥0.10 m/s的分布百分比分別為86.15%、 94.41 %、97.55 %和98.19 %。因此,當入口流速≥0.20 m/s時,97.55%以上的流體流速均大于0.10 m/s,即達到紊流,且隨著入口流速的增加,環(huán)形管內流體流速分布百分比變化減小。
(2) 入口流速對壓力場的影響 環(huán)形固定床環(huán)間距10 mm,管長1.5 m,內部填充玻璃珠直徑為5.86 mm時,不同入口流速在環(huán)形固定床內的壓力變化如圖9所示,L=0.25 m斷面壓差隨流速的變化如圖10所示。從圖9可以看出隨著管長的增加,壓力呈線性減小的趨勢,到達出口處壓力為0;入口流速越大,壓差越大。由圖10看出,在同一斷面處,隨著入口流速的增加,壓差呈拋物線趨勢增大。因此,選取入口速度為0.2 m/s,既能達到管內流速紊流狀態(tài)的要求,又能保證管內壓差適中。
圖7 不同入口流速L/2斷面流速云圖Fig.7 Cloud diagram of velocity at different inlet velocity L/2 sections
圖8 不同入口流速的環(huán)形管內流體流速分布Fig.8 Fluid velocity distribution in annular tube with different inlet velocity
研究以入口流速為0.2 m/s為例,考察了玻璃珠直徑對環(huán)形固定床反應器內流場和壓力場的影響。
(1) 玻璃珠直徑對流場的影響 不同直徑玻璃珠在L/2斷面處的流速云圖如圖11所示。由圖11
圖9 不同入口流速在環(huán)形固定床內的壓力變化Fig.9 Pressure changes of different inlet flow rates in the annular fixed bed
圖10 L=0.25 m斷面壓差隨流速的變化Fig.10 Influence of pressure difference in sectionL=0.25 m on with different flow velocity
圖11 不同玻璃珠直徑L/2斷面流速云圖Fig.11 Cloud diagram of velocity of section L/2 of different glass bead diameters
可以看出,環(huán)形管間隙內流體流速較大,玻璃珠表面一定區(qū)域內流速非常小。這是因為流體流動過程中遇到玻璃珠,玻璃珠與流體之間的黏性導致阻力增大,流速減小,致使玻璃珠周圍流速很小。并且直徑為4.88 mm、5.86 mm、6.93 mm、7.91 mm的玻璃珠其流速在0.5 m/s以上的百分比分別為25.51%、 28.85%、7.34%和8.50%,6.93 mm和7.91 mm玻璃珠環(huán)形管間隙較大,流速大于0.5 m/s的百分比較小,間隙內流速分布較為均勻,說明玻璃珠在阻礙流體流動的同時也有促使流體達到紊流的作用。
不同直徑玻璃珠在入口流速為0.2 m/s時,環(huán)形管內部流體的流速分布百分比如圖12所示,4.88 mm和5.86 mm、6.93 mm和7.91 mm玻璃珠的整體流速分布百分比相似,在0.4~0.5 m/s的范圍,百分比均為最大值,即在0.4~0.5 m/s的范圍內,4.88 mm、5.86 mm、6.93 mm、7.91 mm玻璃珠的流速分布百分比均達到最高值,分別為29.36%、24.80%、39.72%和36.34%。在0.1~0.5 m/s的范圍內,4.88 mm、5.86 mm、6.93 mm、7.91 mm玻璃珠的流速分布百分比分別為72.97%、68.70%、87.27%和87.26%。
圖12 不同玻璃珠直徑的環(huán)形管內流體流速分布Fig.12 Fluid velocity distribution in annular tube with different glass bead diameters
對比圖11與圖12發(fā)現(xiàn),在入口速度為0.2 m/s、玻璃珠直徑為6.93 mm和7.91 mm時,87.27%和87.26%的流體既能保證管內紊流狀態(tài)又能保證管內能量消耗較小。
(2) 玻璃珠直徑對壓力場的影響 不同直徑玻璃珠的壓力隨管長的變化如圖13所示,隨著管長的增大,壓力逐漸減小且呈良好的線性關系,即入口處壓力最大,出口處壓力為0,這符合壓力出口的邊界條件。隨著玻璃珠的增大,壓差逐漸減小,即7.91 mm玻璃珠的耗能最小。這主要是因為玻璃珠越大,環(huán)形管內孔隙率越大,阻力越小,壓差越小。
綜合考慮入口流速、玻璃珠直徑對流場和壓力場的影響,在環(huán)形管結構不變時,入口流速0.2 m/s、玻璃珠直徑7.91 mm為最佳模擬工況。
對數(shù)值模擬優(yōu)化的工況作進一步降解效果評價,驗證數(shù)值模擬的正確性。以苯酚為目標降解物,利用圖1(b)實驗裝置,裝填不同直徑玻璃珠進行了不同入口流速的光催化降解實驗。
圖13 不同玻璃珠直徑的壓力隨管長的變化Fig.13 Pressure changes of different glass bead diameters with different pipe length
實驗條件與測試儀器:苯酚初始質量濃度為20 mg/L,反應容積為30 L,實驗光源為主波長254 nm的低壓汞燈,實驗反應開始前先避光循環(huán)1 h、取3個平行樣,然后開啟光源,運行4 h取樣測試;采用總有機碳(TOC,total organic carbon)分析儀(TOC-Lcpn)測試苯酚降解前后TOC的量。
入口流速對苯酚TOC降解量的影響如圖14所示。由圖14可知,不同直徑玻璃珠隨著入口速度的增加,苯酚TOC礦化量先增加后降低,在速度為0.15~0.20 m/s的范圍內達到最大值,環(huán)形固定床內達到紊流狀態(tài),流體與催化劑充分接觸,反應速率最快,降解效果最好。當速度小于0.15 m/s時,紊流不充分,質量交換速率小,苯酚礦化量小;當速度大于0.2 m/s時,管內流速大,苯酚在玻璃珠表面停留時間短,光照時間短且不充分,苯酚礦化量小。
圖14 入口流速對苯酚TOC降解量的影響Fig.14 Influence of inlet flow rate on degradation of phenol TOC
玻璃珠直徑對苯酚礦化量影響顯著,4.88 mm、7.91 mm玻璃珠對苯酚礦化量效果優(yōu)于5.86 mm、6.93 mm玻璃珠的處理效果,原因是玻璃珠表面光催化界面的流態(tài)、光輻射強度和催化劑量的綜合影響作用,其中7.91 mm玻璃珠在入口速度為0.2 m/s時對苯酚礦化量的效果最好,與數(shù)值模擬的最佳工況結果一致。
根據(jù)實驗和理論壓差結果確定玻璃珠直徑和入口流速與內部阻抗值的關系,然后利用CFD軟件對環(huán)形固定床光催化反應器的流場和壓力場進行了數(shù)值模擬,得到了反應器內流場和壓力場的分布情況。在環(huán)形管結構不變的情況下,通過改變流體入口流速和玻璃珠直徑,分析二者對流場和壓力場的影響,并與實驗結果對比分析。研究結果對環(huán)形固定床反應器的設計及應用具有一定的指導意義。
(1) 確定了玻璃珠直徑、入口流速與內部阻抗值的關系,即y=y0+(A/(w×sqrt(pi/2)))×exp(-2×((v-vc)/w)2),保證數(shù)值模擬結果的準確性。
(2) 玻璃珠直徑為5.86 mm時,管內流體流速和壓降均與入口流速呈正相關,且壓降與流速的平方成正比;當入口流速≥0.2 m/s時,環(huán)形管內97.55%以上的流體流速均大于0.1 m/s,且壓差適中。
(3) 入口流速為0.2 m/s,玻璃珠直徑分別為6.93 mm和7.91 mm時,87.27%和87.26%的流體既能保證管內紊流狀態(tài)又能保證管內能量消耗較小。
(4) 光催化降解苯酚實驗研究顯示,入口流速為0.2 m/s、玻璃珠直徑為7.91 mm時苯酚的降解效果最好,說明了CFD軟件模擬環(huán)形固定床反應器的準確性和可行性。
數(shù)值模擬、理論結果與實驗壓差三者結果吻合較好,說明Ergun方程可用于環(huán)形固定床的壓降計算,也說明了數(shù)值模擬的準確性和可行性。