郭飛宏 王澤宇 仲兆平 王肖祎

(1東南大學(xué)能源與環(huán)境學(xué)院, 南京 210096)(2東南大學(xué)建筑設(shè)計(jì)院有限公司, 南京 210096)

基于球元重建與大渦模擬耦合并行算法的數(shù)值模擬

郭飛宏1王澤宇1仲兆平1王肖祎2

(1東南大學(xué)能源與環(huán)境學(xué)院, 南京 210096)(2東南大學(xué)建筑設(shè)計(jì)院有限公司, 南京 210096)

為了更好地模擬柱形顆粒與石英砂的流動(dòng),利用球元重建方法對(duì)柱形顆粒進(jìn)行球元重建.采用大渦模擬中的SGS(sub-grid scale)亞格子模型,對(duì)柱形顆粒造成的湍流進(jìn)行解析.模擬中耦合并行算法,縮短了模型計(jì)算時(shí)間.基于上述方法,模擬了表觀氣速為1.5 m/s的柱形顆粒與石英砂的混合流動(dòng)過程. 通過流動(dòng)瞬時(shí)圖像和壓力脈動(dòng)分析表明,模擬結(jié)果與實(shí)驗(yàn)結(jié)果相似,表明耦合并行算法能較好地模擬流化床中柱形顆粒與石英砂的流動(dòng).當(dāng)柱形顆粒尺寸較大時(shí),加劇了湍流效應(yīng)對(duì)混合流動(dòng)的影響.模擬的z方向氣體平均速度呈現(xiàn)中間高兩邊低的對(duì)稱分布狀態(tài),這與床內(nèi)的氣泡和柱形顆粒分布有關(guān),符合實(shí)際流化規(guī)律.

柱形顆粒;石英砂;模擬;球元重建;大渦模擬;并行算法

傳統(tǒng)的離散元和計(jì)算流體力學(xué)方法(DEM-CFD)只能對(duì)球形顆粒的流動(dòng)進(jìn)行模擬,然而在實(shí)際的硫化床氣固反應(yīng)過程中,有很多柱形顆粒[1-3].柱形顆粒形狀改變了顆粒與湍流的相互作用,不同形狀顆粒在流場(chǎng)中的碰撞和受力有顯著的差異[4-5].在研究方法上,缺少對(duì)柱形顆粒系統(tǒng)的實(shí)驗(yàn)研究,針對(duì)球形顆粒的研究方法應(yīng)用于柱形顆粒系統(tǒng)時(shí),容易產(chǎn)生很大的誤差.對(duì)流態(tài)化下柱形顆粒的運(yùn)動(dòng)開展數(shù)值模擬也具有一定的難度,如何構(gòu)建柱形顆粒,采用較為統(tǒng)一的數(shù)理模型描述柱形顆粒的受力和運(yùn)動(dòng),是至今仍未完全解決的難題.國(guó)內(nèi)外很多研究者對(duì)DEM-CFD方法進(jìn)行了改進(jìn),Hilton等[6]采用超曲面法建立非球形顆粒模型,模擬了非球形顆粒在流化床中的流動(dòng).Zhong等[7]對(duì)氣固流化床中的柱形顆粒流動(dòng)進(jìn)行了模擬,將柱形顆粒抽象為球元的組合體.Ren等[8-9]對(duì)流化床中玉米狀顆粒流動(dòng)進(jìn)行了模擬,與Zhong等[7]的方法一樣,仍然是將非球形顆粒抽象為球元的結(jié)合體,顆粒所受氣相曳力近似等于球元所受曳力之和.盧洲等[10]采用DEM-CFD方法模擬了彎管中輸送過程的柱狀顆粒.Oschmann等[11]采用球元聚集法和多面體法構(gòu)造柱形顆粒并描述顆粒之間的運(yùn)動(dòng).當(dāng)顆粒尺寸較大時(shí),無法忽視流化床湍流對(duì)氣固流動(dòng)的影響.大渦數(shù)值模擬對(duì)于尺寸大于網(wǎng)格的湍流,采用直接數(shù)值求解[12-14];對(duì)于尺寸小于網(wǎng)格的小尺度的湍流脈動(dòng),采用SGS模型求解.由于大渦模擬能在較短的計(jì)算時(shí)間內(nèi)獲得比雷諾平均模擬更多更準(zhǔn)確的湍流信息.因此,大渦模擬在湍流解析上得到愈來愈廣泛的運(yùn)用.

本文采用球元重建方法,將柱形顆粒表面所受的曳力轉(zhuǎn)化為虛擬球元所受的曳力.采用大渦模型,并引入并行算法,從而實(shí)現(xiàn)柱形顆粒與石英砂流動(dòng)的DEM-CFD并行模擬.

1 基本理論

1.1 球元重建

計(jì)算CFD氣固曳力時(shí),采用直徑為0.8 mm虛擬球元重建一個(gè)長(zhǎng)10 mm、直徑10 mm的柱形顆粒模擬真實(shí)柱形顆粒(見圖1).本文采用0.8 mm的虛擬球元,一方面虛擬球元的直徑與石英砂顆粒直徑大小相同,利于計(jì)算;另一方面,球元之間成正四面體緊密堆積,此時(shí)網(wǎng)格的空隙率約為0.2.當(dāng)CFD計(jì)算網(wǎng)格完全處于柱形顆粒內(nèi)部時(shí),既能保證計(jì)算過程有很好的收斂性,又能借助虛擬球元實(shí)現(xiàn)在柱形顆粒內(nèi)部對(duì)流體形成較大的阻力,使得柱形顆粒內(nèi)部流體速度接近于零.在計(jì)算柱形顆粒與其他顆粒及壁面碰撞時(shí),由直徑為1.5 mm的虛擬球元重建形成一個(gè)中空的柱形殼體,如圖1(c)所示.選用1.5 mm的虛擬球元有利于計(jì)算,同時(shí)保證了柱形顆粒重構(gòu)的準(zhǔn)確性.采用殼體布置方式的原因是:DEM計(jì)算過程中,計(jì)算顆粒碰撞力所用的時(shí)間很長(zhǎng),且計(jì)算用時(shí)隨著顆粒數(shù)量的增加呈指數(shù)級(jí)增加.研究發(fā)現(xiàn),在雙組分DEM-CFD模擬中,可以忽略柱形顆粒內(nèi)部球元的受力狀況,這樣可以縮短計(jì)算時(shí)間.為了盡可能地接近柱形顆粒的表面,球元與球元之間存在重疊.

(a) 柱形顆粒

(b) CFD球元重建

(c) DEM球元重建

為了驗(yàn)證球元重建的有效性,對(duì)柱形顆粒曳力系數(shù)Cd進(jìn)行定量分析.通過改變虛擬球元大小,進(jìn)而改變柱形顆粒內(nèi)部的空隙率,模擬柱形顆粒曳力系數(shù)隨雷諾數(shù)的變化.如圖2所示,柱形顆粒曳力系數(shù)Cd隨著雷諾數(shù)的增加而逐漸減小,并且最后趨于穩(wěn)定;柱形顆粒曳力系數(shù)Cd隨著虛擬球元所占體積分?jǐn)?shù)a的增加而增加.當(dāng)Re<3 000,α=0.7時(shí),曳力系數(shù)與實(shí)驗(yàn)曲線吻合良好；而當(dāng)Re>3 000,α=0.7時(shí),曳力系數(shù)與實(shí)驗(yàn)值有輕微偏離,而α=0.5時(shí)，與實(shí)驗(yàn)曲線吻合良好.本次模擬的截面氣速在1.0～2.0 m/s之間,因此對(duì)應(yīng)的雷諾數(shù)一般小于3 000,所以控制虛擬球元大小使α=0.7,可以獲得較好的流動(dòng)模擬結(jié)果.

圖2 柱形顆粒曳力系數(shù)Cd隨雷諾數(shù)Re的變化

1.2 大渦模擬

對(duì)于湍流的解析,采用大渦模型中普遍采用的SGS(sub-grid scale)亞格子模型[15].在對(duì)傳統(tǒng)N-S方程進(jìn)行過濾后,大渦模型的連續(xù)性方程和動(dòng)量方程為

(1)

-ε

(2)

(3)

(4)

(5)

式中,Cs為SGS常數(shù),在本次模擬研究中取0.13;Δ為特征長(zhǎng)度,根據(jù)孔隙率的不同,其表達(dá)式為

(6)

1.3 氣固作用力

(7)

(8)

本研究中,根據(jù)空隙率的不同分別采用Ergun和en&Yu模型來表示氣固曳力系數(shù)β,表達(dá)式為

(9)

式中,Cd為流場(chǎng)內(nèi)的顆粒曳力系數(shù),由Shiller和Naumann方程以及牛頓曳力法則計(jì)算,即

(10)

式中, 雷諾數(shù)Re的表達(dá)式為

(11)

(12)

式中,αp和αc分別為單個(gè)計(jì)算網(wǎng)格內(nèi)石英砂和生物質(zhì)的體積分?jǐn)?shù);αfic為生物質(zhì)內(nèi)虛擬顆粒的體積分?jǐn)?shù);dp和dfic為石英砂顆粒和虛擬顆粒的直徑.

1.4 離散元方法

Tsuji等[5]提出了DEM-CFD方法,將DEM模型引入到稠密氣固兩相流動(dòng)的模擬研究中.其基本思想是假設(shè)顆粒在碰撞過程中可以有輕微的重疊現(xiàn)象并可以持續(xù)一定的時(shí)間.如圖3所示,采用彈簧、阻尼器和滑移器來表示顆粒與顆粒間或顆粒與壁面間的碰撞.

(a) 法線方向 (b) 切線方向

1—彈簧；2—滑移器；3—阻尼器

圖3 DEM顆粒碰撞示意圖

顆粒發(fā)生碰撞時(shí)其作用力在碰撞點(diǎn)上,這個(gè)力稱為顆粒接觸力Fc.顆粒接觸力被分為法向力Fcn和切向力Fct,法向力使顆粒發(fā)生平動(dòng),切向力由于作用方向線與指點(diǎn)不在同一直線上,會(huì)產(chǎn)生力矩使顆粒發(fā)生轉(zhuǎn)動(dòng).法向力Fcn和切向力Fct由下式計(jì)算：

Fcn=-kndn-ηnvn

(13)

vn=(vr·n)n

(14)

Fct=-ktdt-ηtvt

(15)

vt=vr-(vr·n)n+rwp×n

(16)

式中,dn和dt分別表示法向和切向的相對(duì)位移;vr為顆粒間相對(duì)速度;vt表示碰撞點(diǎn)的相對(duì)滑移速度;kn和kt分別表示法向和切向方向彈簧的彈性系數(shù);ηn和ηt分別為法向和切向的阻尼系數(shù);n為單位法向量;wp為顆粒間相對(duì)角速度.阻尼系數(shù)η可由恢復(fù)系數(shù)e確定,即

(17)

(18)

α=-(1/π)lne

(19)

式中，m為質(zhì)量；k為彈性系數(shù).在計(jì)算柱形顆粒接觸力時(shí),只需計(jì)算柱形顆粒表面的球元與其他顆?；虮诿娴呐鲎?柱形顆粒表面球元碰撞接觸力計(jì)算方法與球形顆粒是相同的.在柱形顆粒與其他顆粒碰撞時(shí),可能只有一個(gè)球元受力,也可能有多個(gè)球元同時(shí)受力.球形顆粒與柱形顆粒的碰撞視為球形顆粒與球元的碰撞.柱形顆?？偟氖芰橹晤w粒內(nèi)所有球元的受力.柱形顆粒的運(yùn)動(dòng)可以認(rèn)為是由球元的平動(dòng)和轉(zhuǎn)動(dòng)疊加產(chǎn)生的.根據(jù)牛頓第二定律,球形顆粒i的平動(dòng)和轉(zhuǎn)動(dòng)方程為

(20)

(21)

式中,mp為顆粒質(zhì)量;vpi為顆粒i的速度矢量;fpi為顆粒i所受氣相曳力;∑Fci為顆粒i所受接觸力總和;Ip為顆粒的轉(zhuǎn)動(dòng)慣量;ωpi為顆粒i的角速度矢量;∑Mpi為顆粒i所受力矩總和.

1.5 并行算法

對(duì)柱形顆粒與石英砂流動(dòng)進(jìn)行模擬,采用了區(qū)域分解技術(shù)和顆粒分配技術(shù)相結(jié)合的方法[16-17].基于信息傳遞的并行編程技術(shù)(message passing interface,MPI)對(duì)CFD求解區(qū)域和顆粒碰撞計(jì)算進(jìn)行多線程并行算法優(yōu)化.在CFD計(jì)算中,連續(xù)相的速度、壓力等重要參數(shù)需要耦合迭代計(jì)算,一個(gè)時(shí)間步長(zhǎng)需要進(jìn)行上千次的迭代計(jì)算,因此連續(xù)相的速度、壓力的計(jì)算是一個(gè)十分復(fù)雜且耗時(shí)的過程,對(duì)整個(gè)計(jì)算區(qū)域按照x方向進(jìn)行區(qū)域劃分.程序運(yùn)行后,在同一時(shí)刻存在N個(gè)線程對(duì)不同計(jì)算區(qū)域進(jìn)行CFD求解,這樣單位時(shí)間內(nèi)進(jìn)行的迭代次數(shù)就會(huì)明顯增加.而每個(gè)計(jì)算域都比原先劃分多了左右2列虛擬邊界作為數(shù)據(jù)緩沖區(qū),這些區(qū)域存儲(chǔ)了與當(dāng)前計(jì)算域相鄰的邊界信息.在一個(gè)迭代計(jì)算步驟中,分區(qū)的邊界信息通過通信句柄MPI_send將計(jì)算域的邊界信息傳輸至虛擬邊界,通過通信句柄MPI_recv從相鄰區(qū)域的緩沖區(qū)獲取邊界信息以保證整個(gè)流場(chǎng)迭代計(jì)算的連續(xù)性.在DEM計(jì)算中,由于顆粒數(shù)量較多,通常都有幾十萬個(gè),因此顆粒碰撞力占用了絕大部分的DEM計(jì)算時(shí)間,需要對(duì)這部分計(jì)算進(jìn)行并行算法優(yōu)化.系統(tǒng)統(tǒng)一對(duì)柱形顆粒中的虛擬球元和石英砂顆粒分配ID號(hào),然后按照ID號(hào)的大小將不同顆粒的計(jì)算任務(wù)分配給不同的線程處理,從而提高計(jì)算效率.

2 模型的建立

2.1 基本流程圖

本文模擬的計(jì)算流程如圖4所示,采用Fortran語言編寫程序.在一個(gè)時(shí)間步長(zhǎng)內(nèi),模擬程序主要進(jìn)行了兩大部分計(jì)算：① 連續(xù)相求解(CFD).基于正交結(jié)構(gòu)化交錯(cuò)網(wǎng)格和有限體積離散方法求解連續(xù)相的流場(chǎng).將動(dòng)量方程離散后采用SIMPLE算法進(jìn)行計(jì)算.首先假設(shè)一個(gè)壓力場(chǎng)已知,求解流場(chǎng)的速度,然后修正壓力和速度,最后進(jìn)行迭代計(jì)算判斷收斂性.② 離散相求解(DEM).首先檢索網(wǎng)格內(nèi)的顆粒,判斷顆粒是柱形顆粒中的球元顆粒還是石英砂顆粒.采用顆粒碰撞軟球模型和球元重建理論,計(jì)算顆粒間的碰撞力,得到柱形顆粒(虛擬球元)和石英砂的運(yùn)動(dòng)軌跡.

圖4 DEM-CFD模擬計(jì)算流程

2.2 計(jì)算區(qū)域及定解條件

構(gòu)建80 mm×40 mm×500 mm的計(jì)算區(qū)域.如果網(wǎng)格尺寸過大,會(huì)使模擬計(jì)算精度大大降低;如果網(wǎng)格尺寸過小,則氣固曳力公式不適用.在CFD計(jì)算中,網(wǎng)格尺寸一般為顆粒尺寸的2.5～4.0倍.在DEM計(jì)算中,網(wǎng)格尺寸的要求相對(duì)較低,一般為顆粒尺寸的1.0～1.5倍.模擬的計(jì)算區(qū)域與實(shí)驗(yàn)中的流化床幾何尺寸一致,在實(shí)驗(yàn)中觀察到流化床上半部分幾乎沒有顆粒的存在,為節(jié)省模擬時(shí)間,在模擬中只計(jì)算床內(nèi)高度小于500 mm的流動(dòng)情況.流化床的幾何尺寸見圖5.流化床模擬計(jì)算中,壁面采用無滑移邊界條件.入口邊界采用速度入口邊界條件,布風(fēng)板上氣孔內(nèi)設(shè)置為速度入口邊界條件,其余部分為壁面的無滑移邊界條件.出口邊界采用常用的壓力出口邊界條件.對(duì)于冷態(tài)模擬,床內(nèi)壓強(qiáng)和氣相速度的初始條件對(duì)模擬結(jié)果的影響很小,本研究設(shè)置初始?jí)簭?qiáng)為0,初始?xì)庀嗨俣葹?.

(a) 主視圖

(b) 側(cè)視圖

2.3 基本參數(shù)

流化床柱形顆粒-石英砂雙組分流動(dòng)模擬主要參數(shù)見表1.

表1 柱形顆粒-石英砂雙組分流動(dòng)模擬參數(shù)

3 結(jié)果與討論

3.1 流動(dòng)圖像分析

模擬結(jié)果采用Tecplot360圖像后處理軟件輸出.圖6為柱形顆粒質(zhì)量分?jǐn)?shù)為10%、表觀氣速為1.5 m/s時(shí),柱形顆粒與石英砂混合流動(dòng)模擬圖.圖7為該工況下的流動(dòng)瞬時(shí)圖.對(duì)比模擬圖和流動(dòng)瞬時(shí)圖可以發(fā)現(xiàn),模擬結(jié)果與實(shí)驗(yàn)結(jié)果得到的圖像較為相似.當(dāng)t=0 s時(shí),柱形顆粒分布在石英砂上面,流化床處于靜止?fàn)顟B(tài),床高約為60 mm.隨著時(shí)間的增加,顆粒開始流化:當(dāng)t=0.3 s時(shí),床層下部形成比較明顯的大氣泡,柱形顆粒和石英砂在氣泡的作用下被吹起,床層“膨脹”的高度為150 mm左右,約為靜止床高的2.5倍,此時(shí)柱形顆粒主要集中在石英砂上部,雖然被氣泡帶到一定高度,但與石英砂并未較好地混合.當(dāng)t=0.5 s時(shí),氣泡逐漸增大和上升,在床層中部形成較大的氣泡區(qū)間,此時(shí)柱形顆粒已經(jīng)到達(dá)一定的高度,最大位置變化不大,但是由于中間氣泡的上升,部分石英砂和柱形顆粒受到中間力的作用,開始向床兩側(cè)的壁面移動(dòng).當(dāng)t=1.2 s時(shí),氣泡破裂成較小氣泡,柱形顆粒和石英砂回落,并再次與上升的氣泡相聚,造成氣泡重復(fù)破裂(見圖6),整個(gè)過程有利于柱形顆粒和石英砂的混合.在0～1.2 s的過程中,氣泡的形成—生長(zhǎng)—破裂、柱形顆粒的上升和下落、柱形顆粒與石英砂的混合程度在逐漸增加.且圖6的氣固流動(dòng)行為與圖7中實(shí)驗(yàn)拍攝的圖像行為十分相似,模擬效果較好.

(a) t=0 s

(b) t=0.3 s

(c) t=0.5 s

(d) t=1.2 s

(a) t=0 s

(b) t=0.3 s

(c) t=0.5 s

(d) t=1.2 s

當(dāng)t=1.0 s時(shí),柱形顆粒與石英砂的流動(dòng)與模擬如圖8所示.此時(shí),較大柱形顆粒的回落會(huì)對(duì)區(qū)域A內(nèi)的大氣泡造成破裂,形成較小氣泡.回落的顆粒再次與上升的氣泡相聚,造成氣泡重復(fù)破裂,促進(jìn)柱形顆粒和石英砂的混合.圖8中的流動(dòng)圖像和模擬圖像均較好地反映出了這一現(xiàn)象.同時(shí),由于柱形顆粒體積較大,四周的石英砂與其碰撞后容易向外側(cè)運(yùn)動(dòng),在柱形顆粒周圍形成較大空隙.這樣,在計(jì)算氣固曳力時(shí)就需要考慮湍流對(duì)氣固流動(dòng)的影響.

(a) 模擬圖 (b) 流動(dòng)圖

3.2 壓力脈動(dòng)分析

圖9為不同床層高度h的壓降模擬值與實(shí)驗(yàn)值對(duì)比.由于最小流化速度為0.5 m/s,因此模擬沒有包含固定床工況(表觀氣速v<0.5 m/s).當(dāng)床高h(yuǎn)=50 mm時(shí),實(shí)驗(yàn)值的床層壓降隨著氣速的增加,先上升后下降并趨于穩(wěn)定,流化后的平均壓降ΔP為1.3 kPa.模擬值主要處于流化階段,床層平均壓降約為1.2 kPa.當(dāng)床高h(yuǎn)=85 mm時(shí),壓降

(a) h=50 mm

(b) h=85 mm

隨著氣速的增加,先上升后下降并趨于穩(wěn)定,流化后的平均壓降約為1.7 kPa,模擬的平均壓降約為1.6 kPa.由于在模擬過程中布風(fēng)條件、床料形狀和尺寸等方面做了一系列簡(jiǎn)化假設(shè),因而模擬結(jié)果與實(shí)驗(yàn)結(jié)果存在一定的偏差,相同氣速條件下床層壓差的偏差量在10%以內(nèi).當(dāng)表觀氣速大于臨界流化風(fēng)速時(shí),模擬方法能較為準(zhǔn)確地模擬分析真實(shí)柱形顆粒與石英砂混合流動(dòng)時(shí)的壓力脈動(dòng).

3.3 湍流分析

在傳統(tǒng)的流化床中,氣速分布均勻,雷諾數(shù)比較小,湍流對(duì)氣固流動(dòng)的影響不很明顯.但是在柱形顆粒和石英砂的混合流動(dòng)過程中,柱形顆粒尺寸比較大,加劇了湍流效應(yīng)對(duì)流動(dòng)的影響,因此有必要對(duì)湍流進(jìn)行分析.大渦模擬的基本思想是對(duì)大尺度湍流用直接數(shù)值模擬方法求解,對(duì)小尺度湍流在添加亞格子尺度應(yīng)力后,采用SGS模型求解.小尺度湍流對(duì)大尺度湍流的影響主要表現(xiàn)在耗散效應(yīng)上,也就是說大渦模擬中的能量耗散是從大渦向小渦串級(jí)傳輸?shù)男问匠霈F(xiàn).SGS湍流黏度的大小在一定程度上反映了小渦對(duì)大渦的影響.利用SGS湍流黏度來定量分析柱形顆粒對(duì)氣體湍流耗散產(chǎn)生的影響.圖10為不同氣速下SGS湍流黏度的分布.圖中紅色為較大的湍流黏度,藍(lán)色為較小的湍流黏度,顏色越深對(duì)應(yīng)的SGS湍流黏度越大.由于噴口附近進(jìn)氣速度較大,雷諾數(shù)較高,在噴口上方很小區(qū)域會(huì)形成湍流黏度最大的區(qū)域.當(dāng)表觀氣速v=1.0 m/s時(shí),氣體速度較低,床內(nèi)柱形顆粒和石英砂的流動(dòng)不很充分.此時(shí),氣體呈氣拴式由床層底部向上緩慢運(yùn)動(dòng),顆粒上下混合不均勻.柱形顆粒聚集在石英砂上部,氣體穿過此區(qū)域后在顆粒層上部形成了較強(qiáng)的尾渦,從而增加了氣體的黏性耗散,SGS湍流黏度較大.當(dāng)表觀氣速v=1.5 m/s時(shí),流化床內(nèi)產(chǎn)生的氣泡變大,混合流動(dòng)加劇,柱形顆粒和石英砂混合較為充分,顆粒在整個(gè)床層內(nèi)分布較為分散,此時(shí)SGS湍流黏度分布相對(duì)于表觀氣速v=1.0 m/s時(shí)較為均勻.

(a) v=1.0 m/s

(b) v=1.5 m/s

3.4 氣相速度

圖11為不同表觀氣速下,不同床層高度z方向氣體平均速度vz的分布.流化床x方向的長(zhǎng)度為80 mm.當(dāng)表觀氣速v=1.0 m/s時(shí),低床層(h=20,40 mm)z方向氣體平均速度受噴口進(jìn)氣影響較大,此時(shí)的氣相速度范圍為1.4～1.8 m/s,平均速度分別為1.6和1.5 m/s.當(dāng)h=80 mm時(shí),z方向氣體平均速度明顯下降,呈現(xiàn)中間低兩邊高的波動(dòng)分布狀態(tài),說明氣體穿過底部床料層后向壁面偏移,并沿著壁面向上運(yùn)動(dòng).這主要是由于柱形顆粒的布置和表觀氣速較低造成的.如圖7(a)所示,柱形顆粒開始集中在x方向中心位置,顆粒本身的重力會(huì)對(duì)底部中心造成附加壓力,這種附加的阻力導(dǎo)致氣體向壓力較低的兩側(cè)壁面運(yùn)動(dòng).當(dāng)表觀氣速較小(1.0 m/s)時(shí),氣體呈氣拴式由床層底部向上緩慢運(yùn)動(dòng),顆粒上下混合不均勻,這種附加壓力會(huì)持續(xù)一段時(shí)間.z方向氣體平均速度的波動(dòng)分布還與布風(fēng)板的設(shè)計(jì)有關(guān).

(a) v=1.0 m/s

(b) v=1.5 m/s

當(dāng)表觀氣速v=1.5 m/s時(shí),床內(nèi)z方向氣體平均速度呈現(xiàn)中間高兩邊低的對(duì)稱分布狀態(tài).當(dāng)床層高度從40 mm增加到200 mm時(shí),z方向氣體最大速度分別從3.0 m/s下降到1.5 m/s;當(dāng)h=40,80 mm時(shí),中間區(qū)域氣體縱向速度是左、右兩側(cè)靠近壁面處的2倍.這是由于在中間床層產(chǎn)生了較大的氣泡(見圖7(b)),床料顆粒受到氣泡的抬升和排擠作用向兩側(cè)回落,因此氣體在中間受到床料的阻力比靠近壁面要小,造成了氣體平均速度呈現(xiàn)中間高兩邊低的對(duì)稱分布狀態(tài).表觀氣速增加后,h=20 mm的低床層同時(shí)受到噴口進(jìn)氣和氣固混合的影響.當(dāng)表觀氣速為1.5 m/s時(shí),z方向氣體平均速度波動(dòng)性明顯加大,靠近壁面x=70～80 mm位置處,氣相速度較低.這主要是因?yàn)楸碛^氣速的增加,加劇了柱形顆粒與石英砂的混合,而這種混合并不均勻.運(yùn)動(dòng)過程中受到多種因素的影響,柱形顆粒和石英砂的分布不均勻,體現(xiàn)在z方向氣體平均速度的波動(dòng)(見圖7(d)).當(dāng)h=200 mm時(shí)氣體已經(jīng)穿過床料層,縱向速度分布較為均勻,接近表觀氣速1.5 m/s.

4 結(jié)論

1) 采用球元重建方法,構(gòu)建了柱形顆粒,利用大渦模擬中的SGS亞格子模型,對(duì)湍流進(jìn)行解析,同時(shí)耦合并行算法,縮短了計(jì)算時(shí)間.

2) 基于球元重建與大渦模擬耦合并行算法對(duì)柱形顆粒和石英砂的混合流動(dòng)進(jìn)行了模擬,模擬結(jié)果與實(shí)驗(yàn)結(jié)果相似.

3) 模擬結(jié)果表明:較大柱形顆粒的上升和回落會(huì)造成氣泡的破裂,促進(jìn)柱形顆粒和石英砂的混合;柱形顆粒尺寸較大,加劇了湍流效應(yīng)對(duì)流動(dòng)的影響,有必要對(duì)湍流進(jìn)行解析;當(dāng)表觀氣速v=1.5 m/s時(shí),床內(nèi)z方向氣體平均速度呈現(xiàn)中間高兩邊低的對(duì)稱分布狀態(tài),這主要與床內(nèi)的氣泡和柱形顆粒分布有關(guān),符合實(shí)際流化規(guī)律.

References)

[1]Khan A A, de Jong W, Jansens P J, et al. Biomass combustion in fluidized bed boilers: Potential problems and remedies[J].FuelProcessingTechnology, 2009, 90(1): 21-50. DOI: 10.1016/j.fuproc.2008.07.012.

[2]Cheng H, Hu Y. Municipal solid waste (MSW) as a renewable source of energy: Current and future practices in China[J].BioresourceTechnology, 2010, 101(11): 3816-3824. DOI:10.1016/j.biortech.2010.01.040.

[3]Hilton J E, Cleary P W. The influence of particle shape on flow modes in pneumatic conveying[J].ChemicalEngineeringScience, 2011, 66(3): 231-240. DOI:10.1016/j.ces.2010.09.034.

[4]Wang C, Zhong Z, Wang X. Microscopic flow characteristics in fluidized bed of cylinder-shaped particles[J].KoreanJournalofChemicalEngineering, 2015, 32(12): 2384-2393. DOI:10.1007/s11814-015-0033-y.

[5]Tsuji T, Higashida K, Okuyama Y, et al. Fictitious particle method: A numerical model for flows including dense solids with large size difference[J].AIChEJournal, 2014, 60(5): 1606-1620. DOI:10.1002/aic.14355.

[6]Hilton J E, Mason L R, Cleary P W. Dynamics of gas-solid fluidised beds with non-spherical particle geometry[J].ChemicalEngineeringScience, 2010, 65(5): 1584-1596. DOI:10.1016/j.ces.2009.10.028.

[7]Zhong W Q, Zhang Y, Jin B S, et al. Discrete element method simulation of cylinder-shaped particle flow in a gas-solid fluidized bed[J].ChemicalEngineering&Technology, 2009, 32(3): 386-391. DOI: 10.1002/ceat.200800516.

[8]Ren B, Zhong W, Chen Y, et al. CFD-DEM simulation of spouting of corn-shaped particles[J].Particuology, 2012, 10(5): 562-572. DOI:10.1016/j.partic.2012.03.011.

[9]Ren B, Zhong W, Jiang X, et al. Numerical simulation of spouting of cylindroid particles in a spouted bed[J].TheCanadianJournalofChemicalEngineering, 2014, 92(5): 928-934. DOI:10.1002/cjce.21900.

[10]盧洲,劉雪東,潘兵.基于CFD-DEM方法的柱狀顆粒在彎管中輸送過程的數(shù)值模擬[J].中國(guó)粉體技術(shù),2011,17(5):65-69.DOI:10.3969/j.issn.1008-5548.2011.05.017. Lu Zhou, Liu Xuedong, Pan Bing. Numerical simulation of the process of cylindrical particles in elbow pipe based on CFD-DEM method [J].ChinaPowderTechnology, 2011, 17(5): 65-69. DOI:10.3969/j.issn.1008-5548.2011.05.017. (in Chinese)

[11]Oschmann T, Hold J, Kruggel-Emden H. Numerical investigation of mixing and orientation of non-spherical particles in a model type fluidized bed[J].PowderTechnology, 2014, 258: 304-323. DOI: 10.1016/j.powtec.2014.03.046.

[12]Yamamoto Y, PotthoffM, Tanaka T, et al. Large-eddy simulation of turbulent gas-particle flow in a vertical channel: Effect of considering inter-particle collisions[J].JournalofFluidMechanics, 2001, 442: 303-334.DOI:10.1017/s0022112001005092.

[13]Koukouvinis P, Gavaises M, Li J, et al. Large eddy simulation of diesel injector including cavitation effects and correlation to erosion damage[J].Fuel, 2016, 175: 26-39. DOI: 10.1016/j.fuel.2016.02.037.

[14]Yin B, Yu S, Jia H, et al. Numerical research of diesel spray and atomization coupled cavitation by large eddy simulation (LES) under high injection pressure[J].InternationalJournalofHeatandFluidFlow, 2016, 59: 1-9. DOI: 10.1016/j.ijheatfluidflow.2016.01.005.

[15]Hoover K A, Foley M G, Heasler P G, et al. Sub-grid-scale characterization of channel lengths for use in catchment modeling[J].WaterResourcesResearch, 1991, 27(11): 2865-2873. DOI: 10.1029/91WR01551.

[16]Li T, Garg R, Galvin J, et al. Open-source MFIX-DEM software for gas-solids flows: Part Ⅱ—Validation studies[J].PowderTechnology, 2012, 220(SI): 138-150. DOI: 10.1016/j.powtec.2011.09.020.

[17]Garg R, Galvin J, Li T, et al. Open-source MFIX-DEM software for gas-solids flows: Part Ⅰ—Verification studies[J].PowderTechnology, 2012, 220(SI): 122-137. DOI: 10.1016/j.powtec.2011.09.019.

Numerical simulation based on sphere reconstruction and large eddy simulation coupled parallel algorithm

Guo Feihong1Wang Zeyu1Zhong Zhaoping1Wang Xiaoyi2

(1School of Energy and Environment, Southeast University, Nanjing 210096, China)(2Architectural Design and Research Institute Co., Ltd., Southeast University, Nanjing 210096, China)

To simulate the mixing flow of cylindrical particles and quartz sand, the spherical element was used to reconstruct cylindrical particles. A sub-grid scale (SGS) model was used to analyze the turbulent flow. The coupled parallel algorithm in simulation reduced the computation time. Based on the above method, the mixing flow process of cylindrical particles and the quartz sand at 1.5 m/s was simulated. The simulation results were in good agreement with the experimental results on instantaneous flow and pressure fluctuation. The result shows that the proposed method can simulate the mixing flow of cylindrical particles and the quartz sand in fluidized bed. When the size of cylindrical particles is larger, the influence of the turbulence on mixing flow is aggravated. The average velocity of gas in thezdirection is simulated, presenting symmetrical distribution of high in the middle and low on both sides. This is related to the distribution of bubbles and cylindrical particles in fluidized bed, thus it is in line with the actual fluidization law.

cylindrical particles; quarts sand; simulation; spherical element reconstruction; large eddy simulation; parallel algorithm

10.3969/j.issn.1001-0505.2017.02.015

2016-07-25. 作者簡(jiǎn)介: 郭飛宏(1986—)，男，博士生；仲兆平(聯(lián)系人)，男，博士，教授，博士生導(dǎo)師，zzhong@seu.edu.cn.

國(guó)家自然科學(xué)基金資助項(xiàng)目(U13161115,51276040)、東南大學(xué)優(yōu)秀博士學(xué)位論文培育基金資助項(xiàng)目(YBJJ1644).

郭飛宏,王澤宇,仲兆平,等.基于球元重建與大渦模擬耦合并行算法的數(shù)值模擬[J].東南大學(xué)學(xué)報(bào)(自然科學(xué)版),2017,47(2):283-290.

10.3969/j.issn.1001-0505.2017.02.015.

TQ051

A

1001-0505(2017)02-0283-08