不同氣速下噴動床噴泉區顆粒特性的模擬分析

張彩虹，黃立新，謝普軍，游 鳳，張耀雷

（中國林業科學研究院林產化學工業研究所，生物質化學利用國家工程實驗室，國家林業局林產化學工程重點開放性實驗室，江蘇省生物質能源與材料重點實驗室，江蘇 南京210042）

噴動床自20 世紀50年代問世以來，已被廣泛應用于燃料燃燒、氣化、物料顆粒干燥及催化聚合等工業過程中。由于在噴動過程中顆粒易于形成規律的循環流動，較隨機而復雜的流化床有更多的優點，在工業過程中應用日益廣泛。噴動床內氣相與反應物的反應時間和停留時間依賴于床體幾何結構和運行數，噴動床設計應滿足不同應用要求。盡管噴射區、環形區和噴泉區內的氣固兩相流流體動力特性已有許多的實驗研究，然而由于噴動床內氣固兩相流動的復雜性及測量條件的限制，人們對床體內運行參數對噴動床動力學的影響至今并不清楚。

為推廣使用具有良好傳熱和傳質效果的噴動床，國內外都做了大量的研究工作。劉偉民等[1]自行設計制造了噴動床糧食干燥機，以探索噴動床用于糧食干燥的實際效果和所存在問題的解決方法。He 等[2-3]設計了噴動床，對噴動床內固體顆粒在噴動區、噴泉區及環隙內空隙率及速率變化等進行了實驗研究，為噴動床的流體力學研究提供了大量的實驗數據。隨著計算機技術的迅速發展，大量科研工作者熱衷于噴動床的模擬研究[4-6]。Hosseini 等[7]采用歐拉-歐拉模型對噴動床內氣固流速度矢量及氣流速率進行模擬，考察了恢復系數及曳力系數的顯著性，經過參數的優化研究，得出Wen-Yu 的曳力模型在無導噴管的模擬中為最佳模型。何玉榮等[8]采用顆粒彈性模量考慮顆粒與顆粒非彈性碰撞產生的顆粒作用力，考慮顆粒-顆粒滑動接觸作用力，采用半理論-半經驗顆粒相黏度模型，建立噴動床內氣固兩相雙流體模型。模擬計算揭示同一氣速下噴射區、環形區和噴泉區內顆粒體積分數、速度場等的分布，并與實驗結果進行比較，在一定范圍內有很好的吻合性，但是沒有對不同氣速對噴動床的影響進行分析。

盡管科研工作者對噴動床進行了大量研究，由于氣固兩項間的復雜變化，至今都沒有對其相互影響完全知曉，模擬值與實驗值間仍然存在較大差距。所以，無論是對噴動床模型的數值模擬還是實驗研究數據都仍然相當缺乏。經過查閱相關文獻發現，噴泉區作為噴動床最大氣固表面接觸區，對其模擬研究較少，一般是對同一氣速下噴泉作簡單的模擬，沒有對這一區域變化進行詳細分析研究。本文作者以自行設計的噴動床為模型，進行了氣固兩相模擬研究，即以噴泉區為研究對象，對其在不同氣速下噴泉高度、噴動形態、顆粒體積分數、顆粒速度的變化及影響進行細微分析，獲得模擬結果，為噴動床設計提供理論依據。

1 邊界條件

本研究將歐拉多相流模型應用噴動床模擬，在下列條件基礎上進行：①離散相與連續相共同存在，相互滲透，各相具有不同的速度、溫度和體積分數；②認為離散相也有連續的速度、溫度和體積分數分布；③兩相流體間存在符合一定規律的能量、動量及質量傳遞。氣固兩相均滿足質量守恒方程。根據Du 等[6]對不同曳力模型在噴動床模擬過程中的研究結果顯示，Gidaspow 模型與實驗結果有很好的一致性，所以選Gidaspow 模型來模擬氣固相間的作用力。采用軟件Fluent 6.3 進行模擬計算。

噴動床模型采用有限體積法對控制方程進行離散，Simple 算法對速度與壓力進行耦合求解，模擬采用兩維軸對稱。網格采用結構式劃分方法，并對其進行網格獨立性驗證，網格增加對模擬結果影響不大。空氣為基本相，給定入口邊界為速度入口，入口處湍流定義方法使用湍流強度與水力直徑，湍流強度5%。氣體不可壓縮。以顆粒（綠豆）作為第二相，進口速度為零。出口采用壓力出口，壓力梯度取為零。對氣相，壁面采用無滑移邊界條件。壁面處顆粒相采用自由滑移邊界條件。假設顆粒間碰撞為彈性碰撞，認為碰撞后恢復系數為0.9 時就達到彈性碰撞，最大體積分數采用Du 等[5]的研究結果，即0.59。

對于最小噴動氣速可由Mathur（1994）的經驗公式獲得，見式（1）。

根據式（1）計算得出ums為0.317 m/s，但在此氣速下噴動床無法形成噴動形態。當氣速增加到0.42 m/s 時，雖然沒有形成噴動態，但顆粒表面顆粒開始松動，床層略有升高，軸線處有微凸狀。故將此氣速定位本模擬實驗的最小噴動氣速，其與計算得出的最小噴動氣速（0.32 m/s）約1.3 倍的關系。模擬基本數據如表1 所示。

2 結果與分析討論

2.1 網格獨立性驗證

由于網格的劃分方法對噴動床的模擬結果有重要的影響，故需對網格的劃分方法進行研究。為了驗證網格劃分方法的可行性，以經典模型He 模型進行驗證[2]。He 模型模擬數據見表2。

采用網格軸向間距為2.5 mm，共550 格；徑向是38 格，從軸中心以等比例1.01 沿徑向增加的方法，其模擬結果如圖1 所示，在0.268 m 床層高度軸心處模擬值與實驗值誤差最大，為18.2%，其它誤差均小于此值。可見，此網格劃分方法在一定范圍內具有一定的可行性[9]。故本文采用此結構化方法對本噴動床進行研究。由于網格的多少會影響模擬結果，需要進行網格獨立性研究，軸向網格間距分別為5 mm、4 mm、2 mm；徑向分別是44 格、55 格、110 格，從軸中心以等比例1.01 沿徑向增加，網格總數如表3 所示。3 種網格下模擬空隙率結果如圖2 中所示。由圖2 可以看出，隨著網格數的增加，在噴動區空隙率趨于均勻，不考慮在環隙區顆粒的變化，55000 格與13750 格空隙率的相對差值在5%內。可見變化值較小，為減小計算工作量，選擇13750 格為本模擬實驗網格數。

表1 模擬數據

表2 He 模型模擬數據

2.2 不同氣速噴動現象模擬

圖1 He 模型噴動區空隙率（h=0.218m）

圖2 噴動區空隙率（h=0.2m）

表3 網格數量

考察了固體顆粒在不同氣速下的噴動情況，采用粒徑為0.003 m 的顆粒，在不同氣速下模擬噴動。由于噴泉區為主相與第二相的主要接觸區，無論是傳質還是傳熱過程大部分都會在這個區域內完成，噴泉區內各種狀態的變化都會對結果產生影響，故本文作者主要對噴泉區內顆粒的變化特性進行了研究。圖3 為噴動氣速在(1.0～1.8)ums下顆粒的噴動形態。從圖3 中可以看出，噴動床中心為低顆粒體積分數的噴射區，在噴射區和床壁面間形成高顆粒體積分數的環形區，在噴動床頂部形成相對較低顆粒體積分數的噴泉區。在氣速為1.0ums時，0.003 m的顆粒均不能起噴，這是由于氣速不能克服床層物料所形成得壓力降；隨著進口氣速的增加，噴泉高度逐漸升高。這與Wu 等[10]觀點一致，氣體速率越大產生的噴泉高度就會越高。圖4 為u/ums氣速在1.0～1.8 條件下噴泉高度曲線，認為當顆粒體積分數小于5×10-7時為噴泉最高點。由圖4 可知，除去1.0～1.2 段，噴泉高度與氣速成很好的線性關系，線性方程為Hf=0.4279(u/ums)-0.23，R2=0.9997。對u/ums=1.25 時進行模擬對方程進行了驗證，模擬穩定后Hf=0.304。由方程計算所得Hf=0.3049，誤差小于1%，可見方程與模擬結果具有很好的吻合性，在u/ums為1.2～1.8 范圍內，氣速與噴泉高度成很好的線性關系，但確切結果還有待于實驗驗證。

2.3 噴泉區顆粒體積分數的變化

由于是對噴泉區顆粒體積分數的分析，故所指床層高度均為除去靜止床層0.3 m 以上的噴泉區高度（ZF）。圖5 為噴動氣速為1.1ums時不同高度顆粒體積分數變化。由圖5 可以看出，固體顆粒體積分數在軸心區域較高，隨著徑向距離的增加顆粒體積分數逐漸下降；這是由于在曳力作用下，高體積分數顆粒與高速氣流一起由噴動區進入噴泉區，沿軸心區域向上運動，在噴泉區內空間區域變大，氣體在徑向的擴散作用變大，曳力逐漸變小，顆粒由加速逐漸變為減速運動，當速度為零時，到達最高位置，顆粒在重力作用下開始向下回落，同時又受到徑向擴散氣流的曳力，在徑向緩慢運動。形成圍繞軸心區的環形回落區，由于沿徑向空間區域逐漸變大，空隙率變大，所以出現體積分數下降的趨勢。由圖5 中還可以看出，在軸心區域，顆粒體積分數隨高度的增加呈現先減小后增大的趨勢，而且噴泉寬度隨著高度的增加逐漸減小。這是由于顆粒在低床層位置時，高濃顆粒被高速氣流由噴動區帶入噴泉區，隨著高度增加，受曳力作用較小的外圍顆粒不斷回落，所以顆粒體積分數逐漸下降；當接近最高點時，上升及下降顆粒都急劇在軸心范圍，所以顆粒體積分數逐漸增大。當噴動氣速為1.2～1.8ums時，顆粒體積分數變化與此相同。這種變化趨勢與He 等[3]的實驗結果一致。

圖3 不同氣速下顆粒的噴動形態（u/ums）

圖4 在不同氣速下噴泉高度曲線

圖5 噴泉區顆粒體積分數分布（u/ums=1.1）

圖6 為噴泉區床高為0.025 m 時，不同氣速下顆粒體積分數的變化。由圖6 可以看出，顆粒體積分數隨著噴動氣速的增加而減小，但是減小幅度不同。噴動氣速由1.1ums增大到1.2ums時，軸線處顆粒體積分數約減少0.057；由1.2ums增大到1.3ums時，軸線處顆粒體積分數約減少0.033；由1.3ums增大到1.5ums時，顆粒體積分數基本不變。由1.5ums增大到1.6ums時，軸線處顆粒體積分數約減少0.029；由1.5ums增大到1.6ums時，軸線處顆粒體積分數又變化很小。在隨噴動氣速增加時，總體趨勢是隨著噴動氣速的增大軸心區顆粒體積分數明顯減小，減小幅度逐漸下降，這與Szafran 等[11]的研究結果相吻合。但是出現兩段區間，顆粒體積分數隨（ZF=0.025 m）噴動氣速變化不明顯。由于沒有發現相關文獻，具體結果還需要實驗驗證及探究起因。由圖6 中還可以看出，隨著噴動氣速的增加，在噴泉區0.025 m 床層高度處，顆粒在徑向的分布寬度逐漸增加。這是由于噴動速度較大時，噴泉高度增加，顆粒回落時相對低噴泉高度運行時間就會增大，擴散氣體作用的徑向運行時間也會增加，所以會逐漸變寬。這也解釋了圖5 中所示的，在噴泉外圍區隨噴泉區床層高度增加噴泉寬度逐漸減小的原因。

圖6 不同氣速相同床層高度顆粒體積分數（ZF=0.025 m）

2.4 噴泉區速度的變化

圖7 為噴動氣速為1.1ums時，噴泉區不同床層高度顆粒速度徑向分布。圖7 中顯示顆粒速度隨床層高度的增加逐漸減小，這是由于顆粒受到重力作用的結果。噴泉核心區和噴泉外圍區運動方向不同，在噴泉核心區顆粒向上運動，速度為正值；而在噴泉外圍區顆粒向下流動，速度為負值。沿半徑方向可以看到，在噴泉區軸心處為最大值，且沿徑向逐漸減小，顆粒速度出現零時為顆粒上升區與下降區的分界點。在噴泉外圍區顆粒向下流動，由分界點開始，在重力作用下做加速運動，速度逐漸增加，且隨床層高度增加而減小。這些變化趨勢與He 噴動床實驗噴泉區徑向速度分布結果基本相符。當噴動氣速為1.2～1.8ums時，顆粒速度變化與此相同。

由圖8中相同床層高度顆粒速度變化（ZF=0.025 m）可知，同一床層高度上的顆粒隨氣速的增加速度增大。這是顯而易見的，由于隨著噴動氣速增加，氣固間曳力也會逐漸增大。但隨著噴動氣速的增加，顆粒速率上升的幅度是不同，而且基本無規律可循。由圖中可見，在軸心處，噴動氣速由1.1ums增大到1.8ums過程中，顆粒速度在兩端變化幅度較大，中間變化較小，最大幅度發生在1.1～1.2ums間，約為0.39；最小幅度在1.6ums增大到1.7ums時，速度變化甚微。這可能與模擬過程中的升力、虛體質量力等其它作用力有關。

圖7 噴泉區顆粒速率徑向分布（u/ums=1.1）

圖8 噴泉區相同高度不同速率下徑向分布（ZF=0.025 m）

3 結 論

（1）噴泉高度隨氣速的增加而增大，且在u/ums=1.2～1.8 范圍為線性關系，線性方程為Hf=0.4279(u/ums)-0.23，R2=0.9997。

（2）固體顆粒體積分數在軸心區域較高，隨著徑向距離的增加逐漸下降，且隨高度的增加呈現先減小后增大的趨勢，而且噴泉寬度隨著高度的增加逐漸減小；顆粒體積分數受速率的影響較大，隨著噴動速率的增大顆粒體積分數趨于減小，但減小幅度不同。

（3）在噴泉核心區顆粒速度隨床層高度增大而減小。噴泉核心區和噴泉外圍區運動方向不同，在噴泉核心區顆粒向上運動，速度為正值；而在噴泉外圍區顆粒向下流動，速度為負值。同一床層高度上的顆粒速度隨氣速的增加而增大，但增大值不同。

符 號 說 明

Dc—— 床體直徑，m

Di—— 入口直徑，m

dp—— 顆粒直徑，m

e—— 恢復系數

H—— 床體高，m

H0—— 靜止床層高，m

ums—— 最小流化速度，m/s

ρg—— 氣相密度，kg/m3

ρs——固相密度，kg/m

θ——底角錐度，（°）

ωs.max—— 最大體積分數

[1] 劉偉民，陳健，徐圣言，等. 導向管噴動床小麥干燥工藝優化研究[J]. 農業工程學報，2003，19（6）：197-200.

[2] He Y L，Qin S Z，Lim C J，et al.Particle velocity profile and solid flow patterns in spouted beds [J]. The Canadian Journal of Chemical Engineering，1994，72（8）：561-568.

[3] He Y L，Lim C J，Grace J R，et al.Measurements of voidage profiles in spouted beds[J].The Canadian Journal of Chemical Engineering，1994，72（4）：229-234.

[4] Wu Z，Arun S A.CFD modeling of the gas-particle flow behavior in spouted beds[J].Powder Technology，2008，183：260-272.

[5] Du W，Bao X，Xu J，et al.Computational fluid dynamics （CFD）modeling of spouted bed：Influence of frictional stress，maximum packing limit and coefficient of restitution of particles [J]. Chemical Engineering Science，2006，61：4558-4570.

[6] Du W，Bao X，Xu J，et al. Computational fluid dynamics (CFD)modeling of spouted bed：Assessment of drag coefficient correlations[J].Chemical Engineering Science，2006，61：1401-1420.

[7] Hosseini S H，Zivdar M，Rahimi R.CFD simulation of gas-solid flow in a spouted bed with a non-porous draft tube [J]. Chemical Engineering and Processing，2009，48：1539-1548.

[8] 何玉榮，陸慧林，劉文鐵，等. 噴動床內氣固兩相流體動力行為的數值模擬[J]. 化工學報，2004，55（2）：290-296.

[9] 黃明輝. 噴動床內氣固兩相流體流動數值模擬[D]. 西安：西北大學，2009.

[10] Wu S W M，Lim C J，Epstein N.Hydrodynamics of spouted beds at elevated temperatures[J].Chem.Eng.Comm.，1987，62：251-268.

[11] Szafran R G，Kmiec A. Periodic flucuations of flow and porosity in spouted beds[J].Trans.Porous.Med.，2007，66：187-200.