基于CFD的不同攪拌組合多相流流場研究

姚晨明， 石秀東， 彭晶鑫

(1.江南大學 機械工程學院， 江蘇 無錫 214122；2.江南大學 江蘇省食品先進制造裝備技術重點實驗室， 江蘇 無錫 214122)

在氣液攪拌釜內，流體動力學環境主要由通氣條件、攪拌轉速和湍流組成[1]。通過攪拌輸入機械能，讓整個發酵過程獲得持續穩定的流場、動量、能量及物質傳遞，從而使相關生物化學反應穩定地進行，達到預期的發酵效果[2]。合理的攪拌釜結構設計能夠在發酵過程中為目標產物生產菌種提供一個有利于其生長繁殖的流體動力學環境。在對攪拌釜進行設計的初期，往往都是結合經驗關聯式通過試驗的方式進行，這樣的方式存在研發成本高、研發周期長等問題[3]。因而，CFD數值仿真技術已經成為生物反應器研究中減少研發時間和成本、進行反應器優化及放大設計的重要手段[4-5]。

謝明輝等[6]在不同溶液中通過CFD數值仿真技術與實驗流體力學相結合的方式進行氣液兩相流研究。研究發現氣液傳質性能依賴于槳型結構形式，不同槳型結構產生不同的流場。槳型結構與物性一起決定了氣泡的動力學和傳質性能。曹毅等[7]對實驗室用攪拌釜進行氣液兩相流數值模擬研究，發現單槳葉攪拌釜提供的流場環境難以滿足聚賴氨酸微生物發酵的需要,于是采用增加槳葉數目和槳葉類型對此問題進行了優化。Wang等[8]對不同工況下多相系統的重要水動力參數，如體積平均總氣持率和時間平均局部氣持率、軸向液速等進行了詳細的模擬和分析，進一步預測了氣泡的尺寸分布，揭示了氣相的獨特性質。

課題組以氣液攪拌釜為研究對象，利用CFD技術進行氣液兩相流非穩態求解，探究多種攪拌組合下攪拌釜內流場特性以及能耗情況，為發酵設備中攪拌器的設計及優選提供參考。

1 數值模擬方法

1.1 物理模型

該攪拌釜主要由帶擋板的3層攪拌槳和橢圓封頭的筒體組成。筒體底部裝有環形氣體分布器，氣體分布器的氣體出口向下。考慮到該模型的復雜性以及氣液兩相流模擬時對網格質量的要求較高，對該模型進行了簡化處理，具體見圖1所示。不同的攪拌槳組合，對微生物發酵性能有很大影響，會影響最終產物的產量。同時，多層攪拌時，底槳是氣液分散性能好壞的關鍵[9]。課題組選取6直葉圓盤渦輪槳(RT-6)、6葉半圓式圓盤渦輪槳(CD-6)、下推式45°的4斜直葉槳共3種槳型組成的4種3層攪拌組合進行模擬，攪拌組合見圖2所示,旋轉方向為順時針旋轉，模型具體結構參數見表1。

圖1 攪拌釜簡化模型Figure 1 Simplified model of fermentation tank

圖2 攪拌組合簡化模型Figure 2 Simplified model of mixing combination

1.2 控制方程

課題組使用Eulerian-Eulerian雙流體模型來模擬計算氣液兩相流問題，不考慮溫度及模擬物料的化學反應，因此忽略組分運輸項和源項，簡化的氣液兩相流連續性方程如下：

(1)

表1 模型主要結構參數

式中：φm表示m相的體積分數，即∑φm=1；ρm為該相的密度值，kg·m-3；um代表該相的速度，m·s-1；m代表液相(l)或氣相(g)。

動量方程如下：

-φmp+·(φmτeff,m)+φmρmg+Fex,m+Mm。

(2)

式中：Fex,m為相間作用力，N；Mm表示主相與次相間由于湍流和相對運動等因素造成的動量交換，kg·m·s-1；p為所有氣液相共享的壓力，Pa；φm,ρm,g為氣液相所受重力,N。

τeff,m為第m相的壓力應變張量，其表達式如下：

(3)

式中：λm和μm是m相的剪切和動力黏度,N·s·m-2；I為單位張量。

在氣液兩相流中，曳力是氣液兩相間的主要作用力，對于本研究的氣液兩相流問題采用Kolev等[10]提出的曳力系數模型，在FLUENT中被稱為通用曳力模型(universal drag)。該模型廣泛適用于包括非球形氣泡等各種氣泡流的氣泡變形情況下的氣泡曳力計算，其計算公式如下：

(4)

(5)

式中：CD表示氣液相間曳力系數；Ad為氣泡橫截面積,cm2；ρl表示液相密度值，kg·m-3；Re為雷諾數；ug和ul表示氣相及液相的速度，m·s-1。

在氣液兩相流中由氣泡徑向分布以及液相湍動引起的力叫做湍流擴散力。這種力的作用是使氣泡在徑向上分布得更加均勻，其計算公式如下：

(6)

式中：CTD為湍流擴散系數；CD為氣液兩相間的曳力系數；δt為湍流Schmidt數；υt為亞網格黏度，m2·s-1；φg和φl表示氣相及液相的體積分數。對于本課題氣液兩相流問題采用Burns-et-al模型[11]。

1.3 網格劃分

課題組在模擬過程中采用FLUENT Meshing劃分多面體網格。多面體網格在劃分復雜模型時相較于非結構化自由四面體網格具有更高的精度，而其數量只有四面體網格數量的1/5～1/3[12]。以攪拌組合A為例，取模型上槳葉周圍一點的液相速度，以及總攪拌功率作為指標，考察了3種不同網格數量對液相速度分量以及攪拌功率的影響，具體如表2所示。從表2中可以看出隨著網格數量的增加，該點的液相速度值和攪拌槳的攪拌功率值逐漸趨于平穩，模擬結果幾乎接近，同時考慮到計算機資源和模擬時間等因素，最終選取86.31萬左右的網格數量。不同的攪拌槳組合，劃分的網格數量有所偏差，課題組選取的網格數量基本維持在86萬左右，對最終模擬結果影響不大。

表2 網格無關性驗證結果Table 2 Grid independence verification results

1.4 計算條件

課題組采用ANSYS FLUENT軟件進行數值仿真。氣液兩相流模擬時采用滑移網格模型法(SM)， Eulerain-Eulerain雙流體模型，湍流模型選擇RNGκ-ε模型，進行瞬態求解。求解過程中監測整體氣相體積分數和攪拌軸力矩系數變化，當氣相體積分數以及力矩系數不再發生明顯變化且趨于穩定時，停止求解，此時得到的結果即為最終結果。求解時采用單一氣泡模型，設定氣泡直徑為2 mm[6-7]。Khopkar等[13]研究發現氣液攪拌釜主體區域中虛擬質量力和升力影響較小，因此氣液兩相之間的作用力不考慮升力和虛擬質量力的影響。曳力模型采用universal drag模型，湍流擴散力模型采用Burns-et-al模型。

模擬物料選取空氣和水。具體參數：水的密度為998.3 kg·m-3，黏度為1.00×10-3Pa·s；空氣的密度為1.225 kg·m-3，黏度為1.79×10-5Pa·s。模擬工況為通氣量1.3 m3·h-1，轉速850 r·min-1。

2 結果與討論

2.1 速度場

圖3為相鄰2擋板中間平面上通過數值模擬得到的液相速度云圖(左)與矢量圖(右)。

通過觀察圖3的速度云圖發現組合B和組合D整體區域速度分布不均勻，低速區域范圍較大，特別是上槳葉上方區域整體速度偏小，不利于物質的傳遞與混合，高速區域主要集中在槳葉葉端附近以及槳葉下方；攪拌組合A以及組合C整體區域速度分布相較于組合B和D來說更加均勻，低速區域范圍較小，整體的流動混合效果較佳，物質傳遞效率較好，高速區域主要分布在槳葉葉端周圍。

圖3 液相速度云圖及矢量圖Figure 3 Liquid phase velocity cloud diagram and vector diagram

如圖3的速度矢量圖所示，攪拌組合A和C的右半部分共有4個循環區域。上槳和底槳都是徑向流攪拌槳，其產生的徑向流首先射向壁面，然后向上或者向下流動，形成典型的徑向雙循環流動，與孫東東等[14]描述一致。同時由于中間軸向流攪拌槳的存在，可以看到上槳的下循環區域與底槳的上循環區域之間存在流動交換，這樣有利于循環與循環之間的物質傳遞。攪拌組合B和D的右半區域共有一大一小2個明顯的循環區域，而上槳葉上方區域未見明顯的循環區域。大循環區域主要位于上槳葉與底槳葉之間，其產生是由于徑向流底槳在上槳及中槳形成的軸向循環區域發生了合并從而形成了大循環區域。

圖4(a)～(c)為液相速度在不同高度處的徑向特征線上的分布情況。截取徑向特征線Z1研究攪拌釜上槳和中槳之間區域流動混合情況；截取徑向特征線Z2研究攪拌釜中槳和底槳之間區域流動混合情況；截取徑向特征線Z3研究攪拌釜釜底區域流動混合情況。從圖中可以看出，在3個不同高度Z1=310 mm、Z2=190 mm、Z3=72 mm的徑向特征線上，4種攪拌組合的液相速度的方向分布規律都與其矢量圖一致。在Z1=310 mm處攪拌組合B和D的液相速度在攪拌槳附近區域相比于組合A和C大，然而在靠近壁面處卻相反，說明在上槳位置采用徑向流攪拌槳能夠為附近流體提供很好的徑向流動，避免釜內壁面附近的物質無法參與循環交換。但是在上槳與中槳之間提供的流體軸向循環速度并不能與軸向流攪拌槳相比。在Z2=190 mm處，整體上攪拌組合B和D的液相速度較大，這是因為上槳及中槳都采用了軸向流攪拌槳，使其提供了較大的軸向循環速度。在Z3=72 mm處，以6直葉圓盤渦輪槳為底槳的組合A和B的速度變化規律相似且速度大小相差不大。同樣地，以6半圓葉圓盤渦輪槳為底槳的組合C和D速度變化規律相似且速度大小相差不大。

圖4 液相速度徑向分布圖Figure 4 Radial distribution diagram of liquid phase velocity

從速度場的分析來看，上槳和中槳采用軸向流攪拌槳形成的流場會得到較大的軸向循環速度，但是速度過于集中在槳與槳之間的區域，而上槳采用徑向流攪拌槳后形成的流場整體速度分布較均勻，沒有過于集中的高速區域，更適合于釜內物質的整體循環混合和傳遞。

圖5所示為液相速度在不同高度處特征線上的分布示意圖。

圖5 特征線示意圖Figure 5 Schematic diagram of characteristic line

2.2 氣相體積分數

圖6為氣相體積分數在相鄰2擋板中間平面上的分布情況。從圖6中可以看出：不同攪拌組合所產生的流型不同，氣體分布情況也不同。釜內氣體在每層槳附近及循環區域內部，氣相體積分數明顯高于其他區域。結合圖3液相速度云圖以及矢量圖的分析可以知道導致這種現象的原因是槳葉后方壓力低，在壓差的作用下，氣泡從高壓區向低壓區運動，使得氣泡容易在該區域內聚集，同時攪拌槳所形成的循環渦使得氣體在循環區域內滯留時間延長，增加了該區域的氣相體積分數。而底槳由于采用了圓盤式的結構設計，使得在底槳下方區域聚集了大量氣體，氣相體積分數高。

圖6 氣相體積分數軸向分布Figure 6 Axial distribution of gas volume fraction

圖7為底槳附近橫截面上的氣相體積分數分布云圖。從圖7中可以看到：6直葉圓盤渦輪槳的槳葉后方有氣泡堆積區域，形成了較大的氣穴，而6半圓葉圓盤渦輪槳槳葉后方的氣體堆積區域較少，形成的氣穴并不大[15]，而其攪拌區域周圍氣體分布更均勻一些，平均氣相體積分數也高于6直葉圓盤渦輪槳的攪拌區域。

圖7 橫截面氣相體積分數分布云圖Figure 7 Cross-sectional gas volume fraction distribution cloud map

綜合上述分析可以知道：整體上，攪拌組合A和C的氣體分布比攪拌組合B與D的更均勻一些，有利于好氧微生物的發酵；6半圓葉圓盤渦輪槳相比于6直葉圓盤渦輪槳在徑向上的氣體分散性能更好一些，能使氣體有效地擴散到周圍，更加均勻地分布。因此總體上攪拌組合C的氣液分散效果更好一些。

2.3 通氣攪拌功率

表3給出了各攪拌組合通過數值模擬以及實際試驗得到的通氣攪拌功率。

表3 通氣攪拌功率Table 3 Aeration mixing power

表3中攪拌組合A和C分別由2個徑向流槳和1個軸向流槳組成，而攪拌組合B和D分別由2個軸向流槳和1個徑向流槳組成。從表中可以得到通氣攪拌功率大小順序為組合A>組合C>組合B>組合D。再結合各攪拌組合選用的槳型，又可以知道在相同的工況和結構參數下，徑向流攪拌槳的功率要比軸向流攪拌槳大，6直葉圓盤渦輪槳的功率要比6半圓葉圓盤渦輪槳大。模擬結果與試驗測得的結果基本相符，誤差率在可接受范圍內。同時數值模擬得到的攪拌功率要比實際工作時的攪拌功率偏小，偏小的原因是數值模擬時未考慮流場中存在的固體區域以及其他阻力等因素。

3 結論

課題組基于CFD技術采用滑移網格模型法(SM)研究了氣-液兩相流下不同攪拌組合的非穩態流場特性，得到了以下結論：

1) 在相同工況條件下，不同的攪拌組合形成的流場特性有很大區別。上槳采用徑向流槳的攪拌組合形成的混合流場整體速度分布更均勻，沒有過于集中的高速區域，更適合于釜內物質的整體循環混合和傳遞。同時上槳和底槳采用6半圓葉圓盤渦輪槳的攪拌組合C氣液分散性能更好一些。

2) 總的功耗方面，攪拌組合A>組合C>組合B>組合D。

課題組的研究對發酵設備中攪拌器的設計及優選有重要的參考價值。