基于bolzmann方法的團(tuán)聚物流固交界面處的自由流動(dòng)

基于bolzmann方法的團(tuán)聚物流固交界面處的自由流動(dòng)

1基于amesim的旋轉(zhuǎn)軸系統(tǒng)優(yōu)化體系統(tǒng)存在于自然界，也存在于人們的日常生活、工業(yè)過(guò)程和生態(tài)環(huán)境中，這與提高人類生活水平、發(fā)展國(guó)民經(jīng)濟(jì)密切相關(guān)。流化床就是典型的顆粒流體系統(tǒng),因其具有溫度分布均勻、傳質(zhì)傳熱效果好和易于操作等優(yōu)點(diǎn),廣泛應(yīng)用石油化工等行業(yè)。由于反應(yīng)器內(nèi)顆粒流體動(dòng)力學(xué)特性的復(fù)雜性和實(shí)驗(yàn)手段的局限性,計(jì)算機(jī)模擬已經(jīng)成為研究顆粒流體系統(tǒng)的重要手段。目前,計(jì)算流體力學(xué)(CFD)是研究顆粒流體系統(tǒng)的有效途徑,而表征流固之間相互作用和動(dòng)量交換強(qiáng)度的相間曳力系數(shù)是其中的關(guān)鍵參數(shù)。由于顆粒在流化床反應(yīng)器內(nèi)的分布是非均勻的,宏觀上呈現(xiàn)軸向上稀下濃和徑向內(nèi)稀外濃結(jié)構(gòu),介觀上呈現(xiàn)濃稀分相的團(tuán)聚物或氣泡結(jié)構(gòu),從均相氣固系統(tǒng)或顆粒沉降實(shí)驗(yàn)中總結(jié)的經(jīng)驗(yàn)關(guān)聯(lián)式往往難以適用。針對(duì)以上非均勻現(xiàn)象,Li等采用多尺度方法對(duì)系統(tǒng)進(jìn)行動(dòng)量和質(zhì)量守恒分析,并提出了實(shí)現(xiàn)模型封閉的穩(wěn)定性條件,建立了能量最小多尺度(energy-minimizationmulti-scale,EMMS)模型。至今,EMMS模型已經(jīng)得到了長(zhǎng)足發(fā)展和廣泛應(yīng)用。EMMS模型中引入的流體繞流和穿流團(tuán)聚物的曳力分配原則對(duì)模型的結(jié)果有重要影響,但至今還沒(méi)有完整的獨(dú)立檢驗(yàn),一定程度上影響了EMMS的深化與擴(kuò)展。本文將對(duì)非均勻顆粒群中的曳力分布進(jìn)行直接計(jì)算,分析該分配原則成立的條件,從而促進(jìn)EMMS的進(jìn)一步發(fā)展和應(yīng)用。2在自然條件下的流體固相過(guò)程中,存在大量固相和稀相間的磁場(chǎng)在工程中通用的曳力表達(dá)式很多,常見(jiàn)的有基于流化床實(shí)驗(yàn)的Wen和Yu公式,基于固定床壓降實(shí)驗(yàn)的Ergun方程,以及Gibilaro等對(duì)Ergun方程的改進(jìn)形式,其中,Rep=ρ|U|d/μ為顆粒雷諾數(shù),μ為流體的動(dòng)力學(xué)黏度,Φ為固相份額,U為表觀流體速度,F=Fd/3πμdU為無(wú)量綱曳力,Fd是流體作用在單顆粒上的平均曳力。EMMS將非均勻的氣固系統(tǒng)分解為稀密兩相,稀相由單顆粒均勻懸浮在流體中組成,密相由離散的當(dāng)量直徑為dcl的團(tuán)聚物組成。為了分析稀相與密相之間的相互作用,模型還將密相團(tuán)聚物等價(jià)成直徑為dcl的“顆?！?并認(rèn)為稀相混合物流體繞流該“顆?！?這將產(chǎn)生部分曳力。相關(guān)曳力與曳力系數(shù)的計(jì)算方式如表1,其中涉及8個(gè)變量,分別是:稀相和密相空隙率εf、εc,稀相和密相中表觀流體速度Uf、Uc,稀相和密相中表觀顆粒速度Upf、Upc,密相團(tuán)聚物尺度和體積份額dcl、f。3基于lbm的流固系統(tǒng)模型LBM是1種介觀數(shù)值方法,直接從離散模型出發(fā),應(yīng)用質(zhì)量、動(dòng)量和能量守恒等基本規(guī)律,在分子運(yùn)動(dòng)論和統(tǒng)計(jì)力學(xué)的基礎(chǔ)上構(gòu)架起宏觀與微觀、連續(xù)與離散之間的橋梁,從1種全新的角度詮釋流體運(yùn)動(dòng)的本質(zhì)問(wèn)題。該方法具有算法簡(jiǎn)單、易于并行、能處理復(fù)雜邊界條件、壓力能夠直接求解等優(yōu)勢(shì),在多相流、湍流、滲流等領(lǐng)域得到了廣泛應(yīng)用。利用LBM求得小于顆粒直徑尺度上流場(chǎng)的詳細(xì)信息,進(jìn)而直接求解復(fù)雜顆粒流體系統(tǒng)中流固之間的相互作用力簡(jiǎn)單易行。目前,LBM中常用單一松弛時(shí)間的Bhatnagar-Gross-Krook(BGK)模型和多松弛時(shí)間(MRT)模型,本文采用BGK模型。為了加快計(jì)算速度,本文還采用了通用計(jì)算GPU(generalpurposeGPU,GPGPU)來(lái)快速執(zhí)行該算法,有關(guān)模型與算法的詳細(xì)討論見(jiàn)文獻(xiàn)。4排圓柱的仿真計(jì)算對(duì)于圓柱呈非均勻分布的雙周期無(wú)限大區(qū)域,如圖1所示,只選擇1個(gè)特征單元進(jìn)行計(jì)算。特征單元的無(wú)量綱量寬度Nx=1536δ,高度Ny=2560δ,圓柱直徑d=20δ,其中δ為L(zhǎng)BM計(jì)算中的網(wǎng)格邊長(zhǎng)。建立笛卡爾坐標(biāo)系,平行于計(jì)算單元寬度的方向?yàn)閤軸,平行于高度的方向?yàn)閥軸,原點(diǎn)在計(jì)算單元的左下角。中心團(tuán)聚物的圓心位于計(jì)算單元的中心,另外4個(gè)1/4團(tuán)聚物的圓心分別位于計(jì)算單元的4個(gè)頂點(diǎn)處。在團(tuán)聚物中心處放置第1個(gè)圓柱,然后左右依次放置圓柱并保持圓柱中心間隔均為D,以第1排圓柱為基準(zhǔn)按照三角陣排布其它圓柱從而形成直徑為dcl的圓形團(tuán)聚物,即如圖1所示的虛線圓,并且圓內(nèi)為密相,圓外為稀相(在這里稀相為無(wú)顆粒的純流體)。通過(guò)改變圓柱的間隔可以調(diào)節(jié)團(tuán)聚物(密相)的空隙率εc,改變團(tuán)聚物內(nèi)圓柱的個(gè)數(shù)N可以調(diào)節(jié)計(jì)算單元的固相份額Φ。這里,團(tuán)聚物份額f、團(tuán)聚物空隙率εc、固相份額Φ的計(jì)算方式為:不可壓縮流體在體力g的驅(qū)動(dòng)下流過(guò)計(jì)算區(qū)域,流固邊界采用無(wú)滑移邊界條件,當(dāng)流體受到的體力與圓柱受到的曳力相等時(shí)即達(dá)到定態(tài)。在LBM計(jì)算中,特征計(jì)算單元內(nèi)單圓柱平均曳力系數(shù)CD為其中,,Fk為團(tuán)聚物中第k個(gè)圓柱受到的曳力,可根據(jù)動(dòng)量交換方法直接求得。CD與無(wú)量綱曳力F的關(guān)系為:對(duì)于Wen和Yu公式、Ergun方程或Gibilaro等對(duì)Ergun方程的改進(jìn)形式只需將式(1)、式(2)或式(3)代入式(6)即可求出無(wú)量綱曳力對(duì)應(yīng)的CD;對(duì)于EMMS密相曳力模型有F=Fc+Fi,稀相曳力模型有F=Ff。特別地,在運(yùn)用EMMS曳力模型計(jì)算單顆粒平均曳力系數(shù)時(shí)需要區(qū)分稀相與密相。首先考慮如圖1所示圓柱呈明顯非均勻分布的情況,團(tuán)聚物即為圖1中的密相區(qū),由于圓柱靜止并且稀相區(qū)無(wú)圓柱,從而稀相區(qū)無(wú)需計(jì)算,表1中的部分參數(shù)及表達(dá)式可以得到簡(jiǎn)化,即εf=1,Upc=0,Usc=Uc,Usi=Uf(1-f),而其它參數(shù)不變。從而算得特征計(jì)算單元內(nèi)單圓柱平均曳力系數(shù)為其次,當(dāng)dcl與Nx的取值使得圓柱呈近均勻分布時(shí),可將特征計(jì)算單元全部看作稀相區(qū),同樣地,表1中的部分參數(shù)及表達(dá)式可以得到簡(jiǎn)化,即εf=1-Φ,Upf=0,Usf=Uf,其它參數(shù)不變,從而特征計(jì)算單元內(nèi)單圓柱平均曳力系數(shù)CD為:5固相-固相的自動(dòng)釋放力模型本文首先計(jì)算了N=31個(gè)圓柱構(gòu)成的團(tuán)聚物在20種圓柱間隔r/d=(0.495,1.11,1.585,1.99,2.34,3.725,4.79,5.685,6.475,7.185,7.84,8.455,9.025,9.57,10.085,10.58,11.05,11.505,11.945,12.37)和Rep在0到10之間的曳力系數(shù),這里r=D-d。為了分析Rep-CD關(guān)系,每種構(gòu)型的CD需要在相同Rep的條件下進(jìn)行比較,在LBM計(jì)算中采用的是體力驅(qū)動(dòng),無(wú)法提前預(yù)知達(dá)到穩(wěn)態(tài)之后的Rep,因此需要插值處理,為了減小插值帶來(lái)的誤差,對(duì)同一種構(gòu)型在Rep處于0到10之間計(jì)算了18種體力驅(qū)動(dòng)情況下對(duì)應(yīng)的CD。圖2給出了圓柱間隔r/d=(0.495,12.37),2種情況下LBM計(jì)算結(jié)果與EMMS曳力公式、Wen和Yu公式以及Gilbilaro公式的Rep-CD關(guān)系曲線,2種圓柱間隔對(duì)應(yīng)的團(tuán)聚物空隙率εc=(0.5047,0.9938)。由圖2可以看出,EMMS曳力公式在圓柱密集排布和近均勻排布2種極限情況下都可以準(zhǔn)確地預(yù)測(cè)相間曳力系數(shù),而其它公式只在圓柱近均勻排布時(shí)才能給出合理的曳力系數(shù)。需要指出的是,EMMS曳力系數(shù)曲線在圖2(a)中是根據(jù)式(7)求得,而圖2(b)則是根據(jù)式(8)求得,因此合理地確定系統(tǒng)被視為非均勻分布的條件是運(yùn)用EMMS曳力公式成功預(yù)測(cè)相間曳力的關(guān)鍵。將團(tuán)聚物視為密相時(shí),其曳力來(lái)自兩部分:流體穿過(guò)密相(“穿流”)所產(chǎn)生的微觀作用力和流體繞過(guò)密相(“繞流”)所產(chǎn)生的介觀作用力,而將團(tuán)聚物視為稀相時(shí)曳力只來(lái)自“穿流”所產(chǎn)生的微觀作用力。為此,對(duì)同一種圓柱排布同時(shí)視為均勻與非均勻狀態(tài),并分別通過(guò)EMMS曳力系數(shù)公式(7)和式(8)分別計(jì)算出相同雷諾數(shù)(Rep=1)時(shí)曳力系數(shù)與圓柱間隔之間的變化關(guān)系。如圖3所示,將團(tuán)聚物視為密相時(shí)隨著圓柱間隔的增加曳力系數(shù)(實(shí)線)不斷增加,因?yàn)閳A柱間隔增大提高了介觀曳力在總曳力中的比例,然而實(shí)際情況是隨著圓柱間隔的增加介觀作用力將逐漸消失。與此同時(shí),將團(tuán)聚物視為稀相時(shí)隨著顆粒間隔的增加曳力系數(shù)(虛線)保持不變。因此,可以將實(shí)線與虛線的交點(diǎn)看做應(yīng)用這2種假設(shè)的分界點(diǎn)εx。由圖3(b)可以看出,引入εx后的EMMS曳力公式與LBM計(jì)算結(jié)果吻合良好,而其它曳力公式在顆粒間隔較低時(shí)會(huì)過(guò)高地預(yù)測(cè)曳力系數(shù)。接下來(lái)保持圓柱個(gè)數(shù)N=31不變,研究顆粒雷諾數(shù)Rep對(duì)εx的影響。如圖4所示,當(dāng)固相份額恒定時(shí),分界點(diǎn)εx隨著Rep的增加而呈指數(shù)規(guī)律降低,擬合函數(shù)為εx=0.98686+0.00329exp(-Rep/7.03718)。根據(jù)EMMS理論,當(dāng)介觀作用力可以忽略或消失時(shí),系統(tǒng)從不均勻狀態(tài)向均勻狀態(tài)過(guò)渡。介觀作用力消失意味著以團(tuán)聚物為整體的繞流消失,曳力全部來(lái)自“穿流”產(chǎn)生的微觀作用力。流體與固體壁面碰撞之后動(dòng)量守恒但會(huì)損失能量,隨著圓柱之間的間隔變小即團(tuán)聚物空隙率變小,穿過(guò)團(tuán)聚物所需要的流體能耗增大。雷諾數(shù)增大時(shí)流體的慣性增大,較易穿越密集的圓柱陣列,從而使εx降低。本文還通過(guò)改變圓柱個(gè)數(shù)N=(19,31,55,97,109,121)研究了固相份額對(duì)εx的影響。如圖5所示,Rep恒定時(shí),分界點(diǎn)εx隨固相份額Φ的增大而呈指數(shù)規(guī)律減小,擬合函數(shù)為εx=0.76707+0.23042exp(-Φ/0.14416)。這是因?yàn)?圓柱個(gè)數(shù)越多則團(tuán)聚物在計(jì)算單元中所占面積越大,從而團(tuán)聚物外的流體通道變窄,使得更多的流體穿過(guò)團(tuán)聚物。通過(guò)以上對(duì)Rep-εx和Φ-εx關(guān)系的分析可以看出,分界點(diǎn)εx同時(shí)是Rep和Φ的函數(shù)。為了更好的分析εx受Rep和Φ的影響,進(jìn)一步研究f-Rep-εx關(guān)系。由于Φ=f(1-εx),這相當(dāng)于研究Φ-Rep-εx關(guān)系。包括前面已經(jīng)計(jì)算的算例,對(duì)團(tuán)聚物份額f=(0.1,0.2,0.3,0.4,0.5,0,6,0.7,0.8)、空隙率εc=(0.80,0.85,0.90,0.95,0.96,0.97,0.98,0.99,0.995,0.998)這80種組合在18種體力驅(qū)動(dòng)情況下共1440個(gè)算例進(jìn)行了計(jì)算。結(jié)果表明,在團(tuán)聚物空隙率εc∈[0.80,0.998]的情況下,團(tuán)聚物份額f<0.2時(shí)對(duì)應(yīng)著非均勻區(qū)域,可以利用EMMS的密相曳力公式Eq.(7)求解曳力系數(shù),f>0.6時(shí)對(duì)應(yīng)著均勻區(qū)域,可以利用EMMS的稀相曳力公式Eq.(8)求解曳力系數(shù);在0.2≤f≤0.6時(shí)存在分界點(diǎn)εx,對(duì)應(yīng)的f-Rep-εx關(guān)系如圖6所示,三者之間的關(guān)系可以表達(dá)為:小于分界點(diǎn)εx的部分對(duì)應(yīng)著非均勻區(qū)域,利用EMMS的密相曳力公式Eq.(7)求解曳力系數(shù),大于分界點(diǎn)εx的部分對(duì)應(yīng)著均勻區(qū)域,利用EMMS的稀相曳力公式Eq.(8)求解曳力系數(shù)。6邊界條件假設(shè)本文利用格子Boltzmann方法研究了能量最小多尺度模型中曳力分配方式在不同顆粒雷諾數(shù)與固相份額時(shí)的合理性。采用GPU加速LBM計(jì)算,可以快速給出流體在小于顆粒直徑尺度上的詳細(xì)流動(dòng)細(xì)節(jié),采用直接反彈邊界實(shí)現(xiàn)無(wú)滑移邊界條