崔智文 王 澤 蔣新宇 趙立豪

清華大學(xué)航天航空學(xué)院工程力學(xué)系, 北京 100084

1 引 言

顆粒兩相流常見(jiàn)于工業(yè)生產(chǎn)、自然環(huán)境和生物醫(yī)學(xué)等領(lǐng)域. 20 世紀(jì)以來(lái), 顆粒兩相流相關(guān)的研究主要集中在球形顆粒問(wèn)題, 且現(xiàn)在依然是兩相流領(lǐng)域中最活躍的方向之一. 在這些研究中顆粒的形狀因素通常被忽略, 而簡(jiǎn)化成球形, 并且研究關(guān)注的重點(diǎn)主要是顆粒在流動(dòng)中的聚集和輸運(yùn)現(xiàn)象. 然而, 在實(shí)際流動(dòng)問(wèn)題中顆粒形狀往往不規(guī)則, 顆粒的非球形特性對(duì)顆粒取向與轉(zhuǎn)動(dòng)行為的影響不可忽略, 比如在紙張制造過(guò)程中紙纖維的取向(Lundell et al. 2011a, H?kansson et al. 2014)以及纖維增強(qiáng)材料制備過(guò)程中纖維的排列(Cheung et al. 2009, Soutis 2005)均會(huì)對(duì)最終成品的力學(xué)性能產(chǎn)生重要的影響; 大量的實(shí)驗(yàn)與數(shù)值模擬結(jié)果表明, 纖維以及聚合物對(duì)壁湍流的減阻和調(diào)制作用與其取向行為密切相關(guān)(Paschkewitz et al. 2004, 2005a, 2005b); 不同環(huán)境溫度和濕度下, 形成的冰晶具有不同的形狀, 會(huì)影響冰晶的碰撞融合以及雨雪的形成(Heymsfield 1977, Jucha et al. 2018); 海洋和湖泊中浮游生物往往具有形狀多樣性等(Pedley & Kessler 1992,Saintillan 2018, Qiu et al. 2019, Ardekani et al. 2017a).

近代有關(guān)非球形顆粒兩相流的理論研究最早可以追溯到Jeffery(1922)在Stokes 流中推導(dǎo)的橢球顆粒受剪切作用的力矩方程, 隨后在Batchelor(1970)和 Brenner(1961, 1963b)的努力下, 逐漸完成了Stokes 流動(dòng)中橢球顆粒的阻力及顆粒對(duì)流體產(chǎn)生的額外應(yīng)力(Batchelor 1970, Brenner 1974)的理論表達(dá)式. 由于當(dāng)時(shí)計(jì)算機(jī)計(jì)算能力的限制, 這些研究主要集中在橢球形顆粒在極低雷諾數(shù)流動(dòng)中運(yùn)動(dòng)的理論或?qū)嶒?yàn)研究. 隨著計(jì)算機(jī)硬件水平的飛速提升以及高精度數(shù)值模擬方法的快速發(fā)展, 復(fù)雜流場(chǎng)中大規(guī)模非球形顆粒運(yùn)動(dòng)的高精度數(shù)值模擬得以實(shí)現(xiàn). 研究的顆粒形狀、尺寸以及流場(chǎng)類(lèi)型日益多樣化. 然而, 相較于球形顆粒兩相流的研究, 非球形顆粒兩相流的數(shù)值模擬研究依然處于起步階段, 而非球形顆粒的各向異性的形狀特性也引入了更多的運(yùn)動(dòng)行為 (比如轉(zhuǎn)動(dòng)、取向等)以及顆粒參數(shù)范圍(比如顆粒形狀等). 這使得非球形顆粒兩相流問(wèn)題變得更加復(fù)雜且難以分析, 進(jìn)而給非球形顆粒兩相流的研究帶來(lái)了較大的挑戰(zhàn). 目前, 有關(guān)非球形顆粒兩相流的數(shù)值模擬研究主要集中以下三個(gè)方面: (1) 單個(gè)顆粒行為的研究, 建立對(duì)非球形顆粒在不同流場(chǎng)中的運(yùn)動(dòng)學(xué)和動(dòng)力學(xué)特征以及顆粒與流場(chǎng)之間的相互作用的基本認(rèn)識(shí), 輔助點(diǎn)顆粒模型的建立; (2) 多顆粒行為研究, 以了解顆粒與顆粒、顆粒與流場(chǎng)間的相互作用, 以及建立相關(guān)模型; (3) 研究大規(guī)模顆粒群在復(fù)雜流場(chǎng) (如湍流) 中的運(yùn)動(dòng)行為, 通過(guò)統(tǒng)計(jì)學(xué)方法描述顆粒在流場(chǎng)中的行為特征, 以建立與完善非球形顆粒的歐拉模型.

目前, 非球形顆粒兩相流的數(shù)值模擬主要基于歐拉-拉格朗日的求解框架. 在該框架下, 流體相的Navier-Stokes 方程在歐拉網(wǎng)格上進(jìn)行求解, 而顆粒相則通過(guò)拉格朗日追蹤方法去描述每一個(gè)顆粒的行為. 顆粒在拉格朗日觀(guān)點(diǎn)下的求解能準(zhǔn)確地描述顆粒的平動(dòng)、取向與轉(zhuǎn)動(dòng)行為, 進(jìn)而可以幫助研究者更加直觀(guān)地認(rèn)識(shí)非球形顆粒復(fù)雜的運(yùn)動(dòng)行為與流場(chǎng)的關(guān)系. 基于該框架, 根據(jù)不同的顆粒計(jì)算方法, 顆粒相的模擬可劃分為點(diǎn)顆粒方法與全分辨顆粒方法. 前者用質(zhì)點(diǎn)替代具體形狀的顆粒, 流體與顆粒之間的相互作用均通過(guò)相應(yīng)的理論模型實(shí)現(xiàn). 而后者需要完全分辨顆粒與流體的界面, 流體與顆粒之間的相互作用通過(guò)全分辨數(shù)值模擬直接實(shí)現(xiàn). 近些年來(lái), 國(guó)內(nèi)外學(xué)者使用這兩類(lèi)方法進(jìn)行了大量直接數(shù)值模擬的研究, 使得人們對(duì)非球形顆粒在復(fù)雜流動(dòng)中的行為規(guī)律有了初步認(rèn)識(shí).

本文針對(duì)非球形顆粒兩相流數(shù)值模擬的研究進(jìn)展進(jìn)行綜述. 文章的主體分為三個(gè)部分: 第一部分介紹非球形顆粒兩相流研究的基本概念、數(shù)值方法和理論模型; 第二部分介紹非球形顆粒在簡(jiǎn)單基本流場(chǎng)(均勻來(lái)流、靜止流場(chǎng)和剪切流)和復(fù)雜湍流場(chǎng)(均勻各向同性湍流和壁湍流)中的運(yùn)動(dòng)學(xué)行為; 第三部分針對(duì)湍流中非球形顆粒取向與轉(zhuǎn)動(dòng)行為、非球形顆粒對(duì)湍流的減阻調(diào)制作用的機(jī)理展開(kāi)討論. 最后, 對(duì)未來(lái)值得研究的工作提出一些建議.

2 數(shù)值模擬方法

在本綜述討論的顆粒兩相流中, 顆粒均為離散相. 本節(jié)將主要介紹歐拉-拉格朗日方法中非球形顆粒的描述方法、理論模型以及相關(guān)的數(shù)值模擬手段. 通常, 流體的求解可以采用有限差分方法、有限體積方法、譜方法等常見(jiàn)的數(shù)值方法. 而顆粒的求解根據(jù)建模的方法可以大致分為點(diǎn)顆粒方法與全分辨顆粒方法. 點(diǎn)顆粒方法用質(zhì)點(diǎn)替代具體的顆粒, 流體與顆粒之間的相互作用均通過(guò)相應(yīng)的力學(xué)模型實(shí)現(xiàn). 而全分辨顆粒方法則完全分辨顆粒的形狀, 流體與顆粒之間的相互作用也通過(guò)數(shù)值模擬直接實(shí)現(xiàn). 兩種方法的使用條件主要由研究的主要問(wèn)題確定. 通常點(diǎn)顆粒方法適合顆粒尺寸小于流體的最小特征尺度, 比如湍流中的Kolmogorov 長(zhǎng)度尺度, 此時(shí)顆粒雷諾數(shù)較小, 顆粒尾部不出現(xiàn)明顯的渦脫落, 顆粒的受力以及顆粒對(duì)流場(chǎng)的影響可以通過(guò)點(diǎn)顆粒模型比較準(zhǔn)確地給出. 而全分辨顆粒方法適合顆粒尺寸大于流體最小特征尺寸, 以捕捉點(diǎn)顆粒模型無(wú)法表征的流體慣性特征和對(duì)流場(chǎng)的反饋影響. 點(diǎn)顆粒方法通常使得單顆粒的計(jì)算量較低, 因此,非常適合大量顆粒 (百萬(wàn)量級(jí)) 的計(jì)算, 但點(diǎn)顆粒的模型適用范圍有限, 例如顆粒尺寸足夠小、顆粒雷諾數(shù)須滿(mǎn)足模型要求等. 而全分辨方法通常可以較為精確地模擬單個(gè)顆粒的行為, 但每個(gè)顆粒的計(jì)算成本高, 尤其是對(duì)于尺寸遠(yuǎn)小于流場(chǎng)尺度的顆粒, 難以實(shí)現(xiàn)大規(guī)模顆粒的計(jì)算. 目前這兩種方法在研究中是相輔相成的關(guān)系, 一方面可以通過(guò)全分辨方法提煉點(diǎn)顆粒的模型, 另一方面可利用完善后的點(diǎn)顆粒模型進(jìn)行大規(guī)模顆粒兩相流的模擬.

2.1 非球形顆粒的基本方程

2.1.1 坐標(biāo)系的定義

為了準(zhǔn)確描述非球形顆粒的行為, 一般需要描述顆粒空間位置以及空間的取向(或姿態(tài)).在本文使用三種參考系來(lái)描述單個(gè)顆粒的運(yùn)動(dòng) (如圖1 所示) , 它們分別是:

(1) 慣性參考系 : 用來(lái)描述連續(xù)流場(chǎng)以及顆粒的平動(dòng)行為.

(2) 顆粒參考系 : 以顆粒質(zhì)心位置為原點(diǎn), 軸分別對(duì)齊顆粒的主軸. 顆粒的轉(zhuǎn)動(dòng)行為通常在該參考系下求解.

(3) 隨動(dòng)參考系 : 以顆粒質(zhì)心位置為原點(diǎn),x′′,y′′,z′′軸分別與慣性參考系中的x,y,z軸平行. 與慣性參考系之間相差一個(gè)平動(dòng)位移.

其中, 顆粒參考系與隨動(dòng)參考系之間通過(guò)一個(gè)旋轉(zhuǎn)張量R進(jìn)行轉(zhuǎn)換, 即

圖1描述顆粒運(yùn)動(dòng)的參考系: 慣性參考系 〈 x,y,z〉 , 顆粒參考系 〈 x′,y′,z′〉 以 及隨動(dòng)參考系〈x′′,y′′,z′′〉

同時(shí), 顆粒的取向、轉(zhuǎn)動(dòng)與速度均由矢量進(jìn)行描述. 對(duì)于矢量而言, 其在慣性參考系與隨動(dòng)參考系下的表達(dá)相同. 為了簡(jiǎn)化標(biāo)記, 在不引入歧義的情況下, 通常用慣性參考系的描述方法來(lái)取代隨動(dòng)參考系, 進(jìn)而對(duì)顆粒的運(yùn)動(dòng)行為以及相關(guān)變量進(jìn)行描述. 因此, 對(duì)于慣性參考系描述的矢量p滿(mǎn)足

其中‘(·)′’代表顆粒參考系下的變量, 下同. 對(duì)于慣性參考系中的二階張量H滿(mǎn)足

2.1.2 平動(dòng)與轉(zhuǎn)動(dòng)方程

根據(jù)牛頓第二定律, 顆粒的平動(dòng)與轉(zhuǎn)動(dòng)方程分別表示為

式中,mp為顆粒質(zhì)量,I為顆粒的轉(zhuǎn)動(dòng)慣量,v為顆粒的平動(dòng)速度,為顆粒的轉(zhuǎn)動(dòng)角速度,F與T分別為作用在顆粒上的合外力與合外力矩. 顆粒的轉(zhuǎn)動(dòng)行為一般在顆粒坐標(biāo)系中進(jìn)行求解, 即

如果顆粒的合外力與合外力矩可以準(zhǔn)確得到, 顆粒的運(yùn)動(dòng)行為可以通過(guò)式(1)與式(2)進(jìn)行求解得到.

2.1.3 顆粒取向求解

從顆粒參考系到慣性參考系的旋轉(zhuǎn)矩陣可由歐拉四元數(shù)E=[e0,e1,e2,e3] 來(lái)表示, 即

歐拉四元數(shù)隨時(shí)間推進(jìn)的方程為

式中

2.1.4 顆粒形狀

非球形顆粒通常具有復(fù)雜的外形. 在非球形顆粒的研究中, 主要使用橢球顆粒模型對(duì)非球形顆粒進(jìn)行近似. 主要的原因包括:

(1) 橢球顆粒是最簡(jiǎn)單的非球形顆粒模型. 近一個(gè)世紀(jì)以來(lái), 科學(xué)家們通過(guò)橢球坐標(biāo)系變換以及相關(guān)橢球諧函數(shù)等數(shù)學(xué)方法得到一些橢球顆粒的作用力的解析表達(dá)式. 因此, 橢球顆粒作為非球形顆粒動(dòng)力學(xué)的理論研究的基礎(chǔ)模型, 擁有比較強(qiáng)大的數(shù)學(xué)工具. 盡管如此, 非球形顆粒兩相流的研究道路依然非常困難, 目前得到的理論模型依然有限.

(2) 大量的直接數(shù)值模擬以及實(shí)驗(yàn)研究表明, 橢球顆粒模型能比較好地捕捉到部分非球形顆粒(非橢球形狀)的運(yùn)動(dòng)特性, 比如圓盤(pán)、纖維、圓柱等. 研究橢球顆粒的運(yùn)動(dòng)行為能夠幫助人們認(rèn)識(shí)形狀的各向異性所帶來(lái)的復(fù)雜動(dòng)力學(xué)行為.

在本文中, 橢球顆粒表示為

2.1.5 無(wú)量綱參數(shù)

在非球形顆粒兩相流中常用的無(wú)量綱參數(shù)如下.

(1) 顆粒雷諾數(shù)Rep: 基于顆粒與流體之間相對(duì)滑移的速度與顆粒直徑定義的無(wú)量綱數(shù), 即Rep=2|uf-v|rp/ν. 其中rp為顆粒的參考半徑(通常選取等效半徑或者顆粒的最大半徑).

(3) 平動(dòng)Stokes 數(shù)St: 用來(lái)反映顆粒的平動(dòng)慣性,St=τp/τf.τp與τf分別表示顆粒的平動(dòng)弛豫時(shí)間與流體特征時(shí)間尺度.

(4) 轉(zhuǎn)動(dòng)Stokes 數(shù)Str: 主要在剪切流中使用, 用來(lái)反映顆粒的轉(zhuǎn)動(dòng)慣性,即Str=ρpsrp2/μ.

圖2橢球顆粒模型示意圖.(a)三軸不等顆粒, (b)桿狀顆粒, (c)碟狀顆粒

(5) 密度比D: 表示顆粒與流體密度之比, 即D=ρp/ρf.

需要注意的是在不同的文獻(xiàn)中雷諾數(shù)的定義稍有差別. 在絕大部分問(wèn)題中特征長(zhǎng)度選為顆粒的直徑, 但在有些理論推導(dǎo)的模型中往往選取顆粒的半徑作為特征長(zhǎng)度. 本文已經(jīng)將該差別進(jìn)行了統(tǒng)一, 但為了使形式與原文獻(xiàn)保持一致以方便讀者查閱, 一些地方保留了原文獻(xiàn)中的定義,相關(guān)的差別將會(huì)在需要的時(shí)候提及.

2.2 點(diǎn)顆粒方法

2.2.1 基本思路

點(diǎn)顆粒, 顧名思義, 就是用質(zhì)點(diǎn)替代復(fù)雜形狀顆粒且同時(shí)考慮顆粒取向與轉(zhuǎn)動(dòng)行為. 如果研究的顆粒的尺寸遠(yuǎn)小于流動(dòng)的特征尺度, 點(diǎn)顆粒是一種非常有效的方法(Maxey 2017). 該方法通常需要適當(dāng)?shù)睦碚摶蚪?jīng)驗(yàn)?zāi)Ｐ蛠?lái)描述流體與顆粒之間的相互作用, 顆粒與流體的相互作用主要受到顆粒相的體積分?jǐn)?shù)和質(zhì)量分?jǐn)?shù)的影響(Balachandar & Eaton 2010). 常見(jiàn)的耦合方式如下.

(1) 單向耦合: 當(dāng)體積分?jǐn)?shù)和質(zhì)量分?jǐn)?shù)很小時(shí), 即顆粒稀疏懸浮, 顆粒相對(duì)流體的反饋?zhàn)饔每梢院雎? 只需考慮流體對(duì)顆粒的作用, 即流體作用在顆粒上的力與力矩;

(2) 雙向耦合: 當(dāng)體積分?jǐn)?shù)或質(zhì)量分?jǐn)?shù)較大時(shí), 顆粒相對(duì)流場(chǎng)的反饋?zhàn)饔貌荒芎雎? 此時(shí)需要考慮顆粒對(duì)流體的作用并進(jìn)行耦合, 比如力耦合(Crowe et al. 1977)、力矩耦合(Andersson et al.2012)以及應(yīng)力耦合(Batchelor 1970, Brenner 1974);

(3) 四向耦合: 當(dāng)顆粒相的體積分?jǐn)?shù)進(jìn)一步增大, 顆粒間的相互作用 (如碰撞, 團(tuán)聚和破碎等) 也需要考慮. 此外, 還需要考慮顆粒與顆粒之間的相互作用力(比如靜電力等).

點(diǎn)顆粒求解的基本流程與思路如下.

(1) 獲取顆粒處未擾動(dòng)流場(chǎng). 當(dāng)一步流場(chǎng)更新后, 需要通過(guò)適當(dāng)?shù)牟逯捣椒?比如線(xiàn)性插值、二階拉格朗日插值等)從流場(chǎng)中獲取顆粒位置處的未擾動(dòng)流場(chǎng)的速度與速度梯度. 對(duì)于非球形顆粒而言, 顆粒位置處的速度梯度信息是非常必要的.

(2) 代入合適的顆粒受力模型, 計(jì)算顆粒的動(dòng)力學(xué)方程.

(3) 如果進(jìn)行雙向耦合, 則需要計(jì)算顆粒對(duì)流場(chǎng)的反饋, 并將反饋力、力矩與應(yīng)力耦合到流場(chǎng)中.

(4) 如果進(jìn)行四向耦合, 則需要進(jìn)行顆粒與顆粒的接觸檢測(cè), 并施加相應(yīng)的碰撞模型或者相互作用模型.

非球形顆粒兩相流中常見(jiàn)的流體對(duì)顆粒的作用力以及顆粒對(duì)流體的反饋力模型將在第2.2.2 和2.2.3 小節(jié)進(jìn)行介紹. 因?yàn)榉乔蛐晤w粒之間的相互作用模型及機(jī)制相對(duì)復(fù)雜, 目前相關(guān)研究工作較少, 本文暫不討論顆粒與顆粒之間的相互作用, 而有關(guān)顆粒與顆粒碰撞的模型可以參見(jiàn)第2.3.4 小節(jié).

2.2.2 流體對(duì)顆粒的作用

流體對(duì)顆粒的作用通常體現(xiàn)為力與力矩的作用.

(1)力的作用

式中,Ks為無(wú)量綱Stokes 阻力系數(shù)矩陣,μ為流體的動(dòng)力黏性系數(shù),uf為顆粒處流體速度,為顆粒速度,rp為顆粒的參考半徑(通常采用最大半徑或者等效半徑). 對(duì)于任意形狀的微小橢球顆粒, 其Stokes 阻力系數(shù)矩陣是對(duì)稱(chēng)正定的(Brenner 1963b). 對(duì)于橢球顆粒模型, 在顆粒坐標(biāo)系下, Stokes 阻力系數(shù)矩陣對(duì)角線(xiàn)上的元素非零, 且表示為

對(duì)于具有回轉(zhuǎn)軸的橢球顆粒而言, 令形狀參數(shù)ka=kb=1 ,kc=λ, 其中λ定義為顆粒回轉(zhuǎn)軸長(zhǎng)度與赤道軸長(zhǎng)度之比. 那么與形狀相關(guān)的參數(shù)α,β,γ和χ可以解析地用λ表達(dá).

由于非球形顆粒的阻力系數(shù)矩陣的各向異性, Stokes 的阻力方向通常不與滑移速度共線(xiàn). 與球形顆粒相比, 非球形顆粒受到的阻力不僅依賴(lài)顆粒與流體之間的速度滑移, 而且還取決于顆粒的取向.

同時(shí), 為了反映顆粒趨近流體速度所需要的特征時(shí)間, 通常使用顆粒的弛豫時(shí)間進(jìn)行描述,記為τp. 一般定義τp/τf為顆粒的Stokes 數(shù), 記為St.當(dāng)St趨近于零時(shí), 意味著顆粒幾乎完全跟隨流體運(yùn)動(dòng). 為了評(píng)估非球形顆粒的弛豫時(shí)間, 在湍流問(wèn)題中, 假設(shè)顆粒在取向隨機(jī)時(shí)(Brenner 1963b,Shapiro & Goldenberg 1993), 評(píng)估顆粒運(yùn)動(dòng)對(duì)流場(chǎng)變化的響應(yīng). 此時(shí), 顆粒的弛豫時(shí)間τp可以定義為

(c) Saffman 升力

在顆粒雷諾數(shù)有限小的情況下, 強(qiáng)剪切會(huì)對(duì)運(yùn)動(dòng)的顆粒產(chǎn)生一個(gè)側(cè)向的升力, Saffman(1965, 1968)最早通過(guò)奇異攝動(dòng)理論得到了球形顆粒在簡(jiǎn)單剪切流中的升力模型, 由此稱(chēng)為Saffman 升力. 該模型預(yù)測(cè)的運(yùn)動(dòng)與實(shí)驗(yàn)中觀(guān)察到的結(jié)果符合較好(Taylor 1923), 因此該模型逐漸得到了推廣與應(yīng)用. McLaughlin (1991), Bagchi 和Balachandar (2002), Magnaudet 與他的合作者們(Legendre & Magnaudet 1998, Magnaudet 2003, Magnaudet et al. 2003)對(duì)球形顆粒的Saffman 升力進(jìn)行了改進(jìn)以適應(yīng)更加廣泛的剪切流場(chǎng). 而Harper 和Chang (1968)最早將Saffman 升力拓展到非球形顆粒, 并發(fā)現(xiàn)相關(guān)的升力系數(shù)矩陣與顆粒形狀無(wú)關(guān). 隨后Fan 和Ahmadi (1995),Feng 和Kleinstreuer (2003), Cui Y 等 (2018, 2019, 2020) 在Harper 和Chang (1968)的基礎(chǔ)上繼續(xù)發(fā)展非球形顆粒的剪切誘導(dǎo)升力模型. 其中Cui Y 等(2020) 最近提出了一個(gè)比較高效求解桿狀顆粒剪切升力的近似模型

然而, 目前的剪切誘導(dǎo)的升力模型與全分辨模擬的結(jié)果依然存在比較明顯的差距(Cui Y et al. 2019). 至于該模型式(14)是否適合碟狀顆粒依然有待研究. 目前在絕大部分點(diǎn)顆粒非球形顆粒兩相流的理論研究問(wèn)題中均未考慮剪切升力模型. 升力對(duì)非球形顆粒的動(dòng)力學(xué)行為的影響依然有待進(jìn)一步討論.

(d) 重力與浮力

當(dāng)考慮顆粒沉降等問(wèn)題, 或顆粒Froude 數(shù)小于1, 即Frp=U/gτp ＜1 , 其中U為流體的特征速度,而g為重力加速度, 此時(shí)需要進(jìn)一步考慮重力和浮力的影響(Sardina et al. 2012, Milici et al. 2014)

式中Fg為顆粒受到的重力與浮力的合力, 與mf分別為顆粒質(zhì)量和顆粒排開(kāi)流體的質(zhì)量,g為重力加速度矢量.

(e) 經(jīng)驗(yàn)公式

然而, 在一些實(shí)際工程問(wèn)題中, 顆粒的雷諾數(shù)通常遠(yuǎn)大于1, 此時(shí)顆粒的受到的流體作用力無(wú)法通過(guò)理論的辦法獲得. 為了能夠繼續(xù)采用點(diǎn)顆粒方法, 就需要通過(guò)實(shí)驗(yàn)或者直接數(shù)值模擬的數(shù)據(jù)分析建立相應(yīng)的經(jīng)驗(yàn)?zāi)Ｐ? 而關(guān)于非球形顆粒受流體力和力矩的經(jīng)驗(yàn)?zāi)Ｐ头浅６? 在本文中不進(jìn)行詳細(xì)展開(kāi),具體可見(jiàn)Ouchene et al. (2015, 2016). 本文主要介紹一種基于橢球顆粒的經(jīng)驗(yàn)?zāi)Ｐ?Zastawny et al. 2012). 不同的經(jīng)驗(yàn)?zāi)Ｐ驮谛问缴隙季哂邢嗨菩? 不同點(diǎn)在于系數(shù)的組合形式. 關(guān)于力的模型, 主要被分為阻力FD與升力FL. 在經(jīng)驗(yàn)?zāi)Ｐ椭? 顆粒阻力與滑移速度共線(xiàn), 而升力作用在與滑移速度垂直的方向. 其對(duì)應(yīng)的阻力系數(shù)與升力系數(shù)定義為

在式(19)與式(20)中的擬合系數(shù)ai與bi可參見(jiàn)(Zastawny et al. 2012). 其他相關(guān)的模型形式類(lèi)似, 具體可參見(jiàn)(Mand? & Roséndahl 2010, H?lzer & Sommerfeld 2008).

(2) 力矩作用

(a) 黏性力矩

對(duì)于中心對(duì)稱(chēng)(即具有三個(gè)對(duì)稱(chēng)面)的非球形顆粒, 在Stokes 流的情況下, 其受到的黏性力矩為(Brenner 1974)

上式中,KR為旋轉(zhuǎn)阻力系數(shù)張量,Θ為剪切阻力系數(shù)張量, 其中Θ為三階張量,KR與Θ均與顆粒的形狀有關(guān). 對(duì)于任意三軸不等橢球顆粒而言, 在顆粒參考系中受到的黏性力矩為(Jeffery 1922)

對(duì)于軸對(duì)稱(chēng)橢球顆粒, 形狀有關(guān)的參數(shù)α,β與γ由表1 的表達(dá)式給出.

表1 軸對(duì)稱(chēng)橢球顆粒相關(guān)的形狀參數(shù)

(b) 流體慣性力矩

(c) 經(jīng)驗(yàn)公式

對(duì)于高雷諾數(shù)問(wèn)題, 可將力矩分解成俯仰力矩Tp(pitch torque)與旋轉(zhuǎn)力矩TR(rotation torque). 俯仰力矩主要是當(dāng)顆?；剞D(zhuǎn)軸的方向與流場(chǎng)速度方向不共線(xiàn)時(shí)產(chǎn)生, 而旋轉(zhuǎn)力矩為顆粒相對(duì)于流體旋轉(zhuǎn)時(shí)流場(chǎng)對(duì)顆粒產(chǎn)生的力矩(Zastawny et al. 2012). 在經(jīng)驗(yàn)公式模型中, 通過(guò)定義俯仰力矩系數(shù)為

在Zastawny 等(2012)的工作中, 該系數(shù)與雷諾數(shù)以及攻角的關(guān)系為

上式中擬合系數(shù)c1～c10可參見(jiàn)Zastawny 等(2012)的工作. 旋轉(zhuǎn)力矩系數(shù)定義為

式中Ωs為流體與顆粒之間的相對(duì)滑移角速度, 即Ωs=Ω-ωp. 該系數(shù)與顆粒旋轉(zhuǎn)雷諾數(shù)的關(guān)系為

式中顆粒旋轉(zhuǎn)雷諾數(shù)為ReR=4rp2|Ωs|/ν, 系數(shù)r1～r4參見(jiàn)Zastawny 等(2012)的工作. 其它相關(guān)的擬合形式的力矩形式可參見(jiàn)文獻(xiàn)(Ouchene et al. 2015, 2016; Yin et al. 2003).

(3) 其他受力模型

根據(jù)經(jīng)典的Maxey-Riley 小球受力模型(Maxey & Riley 1983), 非定常流動(dòng)中點(diǎn)顆粒模型還需要進(jìn)一步考慮壓力梯度力、附加質(zhì)量力、Basset 歷史力效應(yīng)項(xiàng)等. 當(dāng)顆粒密度遠(yuǎn)大于流體密度時(shí), 壓力梯度力與附加質(zhì)量力的作用可以忽略不計(jì). 對(duì)于顆粒密度接近流體密度或者小于流體密度時(shí), 壓力梯度力以及附加質(zhì)量力的作用不可忽略. 在當(dāng)前綜述的工作中描述的點(diǎn)顆粒均是密度比遠(yuǎn)大于1 的情況, 暫不考慮壓力梯度力與附加質(zhì)量力的影響. Basset 力反映了顆粒在加速過(guò)程中邊界層發(fā)展滯后所導(dǎo)致的“歷史”效應(yīng)的作用, 但往往因?yàn)槠湓谟?jì)算中實(shí)現(xiàn)困難以及該項(xiàng)僅在具有較大滑移速度變化時(shí)比較重要, Basset 力通常在計(jì)算中被忽略. 近期的研究文章(Olivieri et al. 2014, Daitche 2015, Prasath et al. 2019)表明Basset 力在顆粒輸運(yùn)中的重要影響. 然而,Basset 歷史力效應(yīng)項(xiàng)是基于球形顆粒得到的, 對(duì)于非球形顆粒而言, 目前還沒(méi)有針對(duì)一般橢球顆粒的模型. Lawrence 和Weinbaum (1986)通過(guò)攝動(dòng)理論對(duì)近球形顆粒的小偏心率進(jìn)行展開(kāi), 推導(dǎo)了近球形橢球顆粒的歷史力模型, 并發(fā)現(xiàn)對(duì)于非球形顆粒而言, 除了Basset 歷史效應(yīng)力外還存在新的歷史力項(xiàng). 相關(guān)的推導(dǎo)與解釋可參見(jiàn)論文Lawrence 和Weinbaum (1986, 1988). 目前,在大部分非球形點(diǎn)顆粒模擬中, 依然不考慮歷史力效應(yīng), 而歷史力對(duì)非球形顆粒的影響還有待進(jìn)一步分析.

2.2.3 顆粒對(duì)流體的作用

顆粒對(duì)流體的反饋?zhàn)饔弥饕ㄟ^(guò)力耦合、力矩耦合和應(yīng)力耦合等方法實(shí)現(xiàn), 對(duì)應(yīng)的實(shí)現(xiàn)方法如下.

(1) 力耦合

力耦合是基于牛頓第三定律的傳統(tǒng)雙向耦合方法, 廣泛應(yīng)用于慣性顆粒與流體相互作用的數(shù)值模擬中. 該方法的基本思想是將顆粒所受到的流體作用力作為顆粒對(duì)流體的反饋力. 反饋力的實(shí)現(xiàn)是將某個(gè)網(wǎng)格內(nèi)所有顆粒的反作用力 -F進(jìn)行體積平均 (或者加權(quán)平均) 以得到作用于流體相的歐拉網(wǎng)格點(diǎn)的體積力( Crowe et al. 1977, Maxey et al. 1997, Boivin et al. 1998, Vreman 2015 )

其中,Δ為網(wǎng)格體積,Fil為該網(wǎng)格體積內(nèi)第l個(gè)顆粒所受流體作用力,np為該網(wǎng)格體積內(nèi)的顆粒數(shù). 將體積力fip添加到流體動(dòng)量方程的右端項(xiàng), 便實(shí)現(xiàn)了基于力耦合方法的雙向耦合.

(2) 力矩耦合

Andersson 等 (2012)在力耦合的基礎(chǔ)上提出了力矩耦合方法, 認(rèn)為除了將顆粒的點(diǎn)力作為反饋力外, 還應(yīng)該考慮顆粒所受力矩對(duì)流體的反作用. 和力耦合方法相似, 單個(gè)顆粒受到流體施加的力矩N, 則該顆粒對(duì)周?chē)黧w產(chǎn)生反饋力矩 -N, 其可表示為作用于流體微團(tuán)的表面力. 體積為Δ的網(wǎng)格內(nèi), 顆粒產(chǎn)生的反對(duì)稱(chēng)應(yīng)力為

(3) 應(yīng)力耦合

不同于大密度比的慣性顆粒, 小尺寸的近中性懸浮顆?？梢暈闊o(wú)慣性顆粒, 因此顆粒對(duì)流體的反饋力和力矩可以忽略. 但剛性顆粒由于不可變形, 還會(huì)通過(guò)額外的顆粒應(yīng)力影響附近流體.稀疏懸浮液中的這種應(yīng)力可以通過(guò)對(duì)應(yīng)力子(stresslet)進(jìn)行局部體積平均得到(Batchelor 1970).盡管應(yīng)力子不像力矩和點(diǎn)力那樣與顆粒的動(dòng)力學(xué)直接相關(guān), 但應(yīng)力子實(shí)質(zhì)上代表剛性顆粒對(duì)流體應(yīng)變運(yùn)動(dòng)的抵抗 (Guazzelli & Morris 2011). Brenner (1974)推導(dǎo)了一般軸對(duì)稱(chēng)顆粒的應(yīng)力子形式, 這里只簡(jiǎn)單介紹軸對(duì)稱(chēng)橢球顆粒的對(duì)應(yīng)表達(dá)式

2.2.4 無(wú)慣性橢球顆粒模型

式中,vs為顆粒的Stokes 沉降速度. 對(duì)于小慣性顆粒的沉降問(wèn)題, 顆粒的速度可以通過(guò)顆粒當(dāng)?shù)氐牧鲌?chǎng)速度加上沉降速度獲得(Qiu et al. 2019). 若不考慮重力作用, 顆粒近似跟隨流體軌跡線(xiàn).

此時(shí)顆粒的轉(zhuǎn)動(dòng), 也直接由當(dāng)?shù)氐乃俣忍荻扰c顆粒取向決定

類(lèi)似慣性顆粒的求解, 顆粒的取向相關(guān)的歐拉四元數(shù)由式(5)推進(jìn).

對(duì)于軸對(duì)稱(chēng)的橢球顆粒, 顆粒的取向行為也可以直接通過(guò)求解顆?；剞D(zhuǎn)軸方向方程獲得, 即

2.3 全分辨顆粒方法

對(duì)于尺寸大于流動(dòng)最小尺度 (kolmogorov 尺度) 的顆粒, 以及一些形狀不規(guī)則的顆粒, 它們與流體的相互作用就難以再用點(diǎn)顆粒的方法描述了, 這時(shí)就需要借助全分辨的方法來(lái)對(duì)顆粒兩相流進(jìn)行模擬. 這種模擬方法在顆粒表面施加無(wú)滑移邊界條件, 完全分辨出顆粒的形狀、尺寸,進(jìn)而分辨出顆粒周?chē)牧鲌?chǎng)以及顆粒-流體間的相互作用. 基于完全分辨的方法研究顆粒兩相流能夠最準(zhǔn)確地捕捉顆粒-流體之間的相互作用(Tenneti & Subramaniam 2014), 同時(shí)也可以用于對(duì)基于點(diǎn)顆粒的方法進(jìn)行建模、驗(yàn)證, 具有重要的學(xué)術(shù)價(jià)值和研究意義.

完全分辨的顆粒兩相流問(wèn)題中顆粒通常被視為剛體, 其與流體的相互作用就可以描述為流體中自由運(yùn)動(dòng)剛體的流固耦合問(wèn)題, 如圖3 所示, 其中Si代表顆粒表面, 小箭頭代表顆粒運(yùn)動(dòng)方向. 對(duì)于該系統(tǒng), 流體相用基于歐拉觀(guān)點(diǎn)的不可壓Naiver-Stokes(NS)方程描述, 顆粒相由基于拉格朗日觀(guān)點(diǎn)的動(dòng)力學(xué)方程(1)和(2)描述, 在不考慮其他外力 (如電磁力、電場(chǎng)力等) 的情況下,方程中的右端項(xiàng)合外力、合外力矩的形式為F=FH+Fint,T=TH+Tint, 其中FH,TH表示流體對(duì)顆粒的力和力矩,Fint,Tint表示顆粒之間相互碰撞的力和力矩. 而顆粒表面為流體提供了無(wú)滑移邊界條件

式中uF代表顆粒表面的速度, 和ω分別為顆粒中心點(diǎn)平動(dòng)速度和轉(zhuǎn)動(dòng)角速度,X為顆粒表面點(diǎn)的坐標(biāo), 而Xc為顆粒中心點(diǎn)坐標(biāo). 所以全分辨方式模擬顆粒兩相流自然是雙向 (不考慮顆粒碰撞) 或四向 (考慮顆粒碰撞) 耦合系統(tǒng).

對(duì)于全分辨顆粒兩相流問(wèn)題的求解主要有兩類(lèi)方法, 即貼體網(wǎng)格類(lèi)方法和非貼體網(wǎng)格類(lèi)方法(Haeri & Shrimpton 2012, Kim & Choi 2019). 貼體類(lèi)方法的思想很簡(jiǎn)單, 即在每一時(shí)刻都根據(jù)顆粒的位置劃分流體網(wǎng)格, 并且根據(jù)顆粒的位置和當(dāng)前運(yùn)動(dòng)設(shè)置流場(chǎng)的邊界條件, 然后在流體區(qū)域內(nèi)求解NS 方程. 到下一個(gè)時(shí)間步, 顆粒的位置發(fā)生變化, 則需要對(duì)空間網(wǎng)格進(jìn)行重新劃分以對(duì)流動(dòng)進(jìn)行求解, 如圖4(a)所示. 這類(lèi)方法的典型代表就是任意-拉格朗日-歐拉法 (ALE) .而后一類(lèi)方法則通常需要兩套網(wǎng)格, 即描述流場(chǎng)的背景笛卡爾網(wǎng)格和描述顆粒表面的拉格朗日網(wǎng)格, 如圖4(b)所示的這類(lèi)方法中, 背景的笛卡爾網(wǎng)格并不會(huì)因?yàn)轭w粒的運(yùn)動(dòng)而發(fā)生改變, 因此不需要流體網(wǎng)格隨顆粒運(yùn)動(dòng)的不斷重新生成. 但是需要通過(guò)一些特殊的手段來(lái)處理流體-固體界面的邊界條件. 這類(lèi)方法的典型代表包括浸沒(méi)邊界方法 (immersed boundary method, IBM), 虛擬區(qū)域方法 (fictitious domain method, FDM) 和格子-玻爾茲曼方法 (lattice Boltzmann method,LBM). 貼體類(lèi)方法在追蹤運(yùn)動(dòng)的顆粒的過(guò)程中要不斷更新求解流體方程的網(wǎng)格, 這個(gè)過(guò)程需要很大的計(jì)算量, 所以這類(lèi)方法對(duì)顆粒兩相流的數(shù)值模擬應(yīng)用目前還僅限于二維(Haeri & Shrimpton 2012), 這里不再展開(kāi). 而非貼體類(lèi)方法省去了網(wǎng)格更新過(guò)程的計(jì)算量, 同時(shí)背景的笛卡爾網(wǎng)格的簡(jiǎn)潔性使得流體相NS 方程求解過(guò)程中的一些快速算法得以應(yīng)用. 這兩點(diǎn)使得該方法的計(jì)算效率得到顯著提升, 因而這類(lèi)方法是目前模擬運(yùn)動(dòng)物體流固耦合 (如顆粒兩相流) 問(wèn)題的主流方法. 以下對(duì)非貼體類(lèi)的浸沒(méi)邊界方法和格子-玻爾茲曼方法進(jìn)行展開(kāi)介紹. 由于虛擬邊界法與浸沒(méi)邊界法思路相似, 這里不再展開(kāi)描述, 相關(guān)詳細(xì)過(guò)程可參考(Yu & Shao 2007b, 2010).

2.3.1 浸沒(méi)邊界方法

浸沒(méi)邊界方法最早是由Peskin 提出用來(lái)模擬心臟瓣膜流固耦合問(wèn)題(Peskin 1972, 1977). 如圖5 所示, 浸沒(méi)邊界法使用背景的歐拉網(wǎng)格和物體表面的拉格朗日網(wǎng)格分別描述流體和物體(Griffith & Patankar 2020). 對(duì)于流體-固體界面的邊界條件, 該方法在NS 方程動(dòng)量方程中施加一個(gè)額外的體積力項(xiàng)來(lái)實(shí)現(xiàn), 即將原本NS 方程的動(dòng)量方程修改為

圖3全分辨顆粒兩相流問(wèn)題示意圖

圖4網(wǎng)格方法示意圖. (a)貼體網(wǎng)格類(lèi)方法, (b)非貼體類(lèi)網(wǎng)格方法

圖5浸沒(méi)邊界方法的歐拉網(wǎng)格和拉格朗日網(wǎng)格

式中f稱(chēng)為浸沒(méi)邊界力 (IB 力) . 適當(dāng)?shù)厥┘覫B 力后就可以使得流固無(wú)滑移邊界條件即式(39)得到滿(mǎn)足. 根據(jù)IB 力的施加方法, 浸沒(méi)邊界法可以分為兩類(lèi), 即連續(xù)力法和離散力法(Mittal &Iaccarino 2005). 連續(xù)力法中的IB 力先在物體表面的拉格朗日網(wǎng)格上計(jì)算, 然后再通過(guò)離散的Delta 函數(shù)“傳播”(spread)到背景歐拉網(wǎng)格上. 根據(jù)IB 力的具體計(jì)算方法, 這類(lèi)方法包括Peskin 最早提出浸沒(méi)邊界法時(shí)的經(jīng)典方法(Peskin 2002)、基于彈簧模型的方法(Lai & Peskin 2000)、罰函數(shù)法(Huang et al. 2011)、多孔介質(zhì)模型法(Angot et al. 1999)、反饋力法(Goldstein et al. 1993)、直接力法(Uhlmann 2005)和投影法(Taira & Colonius 2007)等. 離散力法則是對(duì)沒(méi)有IB 力的原始NS 方程先進(jìn)行離散, 然后再根據(jù)流固界面邊界條件在背景歐拉網(wǎng)格上直接施加IB 力. 這類(lèi)方法又包含直接力法(Fadlun et al. 2000)、鬼點(diǎn)法(Majumdar et al. 2001)和切割網(wǎng)格法(Udaykumar et al. 1996)等. 連續(xù)力法中物體邊界的識(shí)別是模糊的, 離散力法中邊界的識(shí)別是清晰的, 但由于離散力法在處理運(yùn)動(dòng)物體時(shí)會(huì)由于邊界位置改變導(dǎo)致插值過(guò)程不連續(xù),進(jìn)而導(dǎo)致虛擬的IB 力振蕩(Mittal & Iaccarino 2005), 所以在處理顆粒兩相流問(wèn)題中, 通常使用連續(xù)力法. 這里以顆粒兩相流全分辨模擬中廣泛使用的直接力法(Uhlmann 2005, Breugem 2012)為例介紹浸沒(méi)邊界法的具體實(shí)施步驟.

首先對(duì)不施加IB 力的NS 方程動(dòng)量方程進(jìn)行離散得到預(yù)測(cè)速度u*

隨后, 對(duì)顆粒的運(yùn)動(dòng) (包括平動(dòng)和轉(zhuǎn)動(dòng)) 利用顆粒動(dòng)力學(xué)方程進(jìn)行推進(jìn)

通常, 對(duì)于非球形顆粒, 轉(zhuǎn)動(dòng)方程在顆粒參考系下求解, 詳見(jiàn)2.1 節(jié).

2.3.2 格子玻爾茲曼方法

格子玻爾茲曼方法 (lattice boltzmann method, LBM) 也是全分辨顆粒兩相流問(wèn)題數(shù)值模擬的有效手段. 格子玻爾茲曼方法本質(zhì)上是一種流體求解器. 但是與傳統(tǒng)的有限差分、有限體積等基于對(duì)描述宏觀(guān)流體運(yùn)動(dòng)的NS 方程進(jìn)行離散的流體求解器不同, 格子玻爾茲曼方法是從介觀(guān)角度的氣體動(dòng)理學(xué)Boltzmann 方程出發(fā)來(lái)描述流體動(dòng)力學(xué)(Chen & Doolen 1998), 即

方程(47)實(shí)際描述的是流體系統(tǒng)中“流體粒子”的運(yùn)動(dòng), 式中的fi是流體粒子的在第i 個(gè)方向的速度分布函數(shù),x代表“晶格” (lattice) 的位置,Ωi代表碰撞算子,i為晶格的方向指標(biāo) (通常對(duì)于二維問(wèn)題可取9 個(gè)方向, 對(duì)于三維問(wèn)題可取15,19 或27 個(gè)方向, 如圖6 所示) , |ei|=dx/dt. 求解方程 (47) 涉及到“對(duì)流”(streaming)和“碰撞”(collision)兩步, 常用Bhatnagar-Gross-Krook(BGK)單弛豫時(shí)間的碰撞算子(Chen et al. 1992, Qian et al. 1992)

式中, 上標(biāo)“eq”代表平衡態(tài)的分布函數(shù), 具體形式可見(jiàn)參考文獻(xiàn)(Chen & Doolen 1998). 得到方程 (47) 的解后, 對(duì)于不可壓流動(dòng)問(wèn)題, 基于Chapman-Enskog 展開(kāi)(Chen & Doolen 1998), 就可以得到流體力學(xué)中的宏觀(guān)量 (如密度、壓強(qiáng)、速度等) . 關(guān)于LBM 的具體求解步驟和推導(dǎo)過(guò)程,可以參見(jiàn)相關(guān)綜述文章(Chen & Doolen 1998, Aidun & Clausen 2010)或相關(guān)書(shū)籍可參見(jiàn) (何雅玲等 2009) .

LBM 方法較之于傳統(tǒng)的方法有以下三個(gè)優(yōu)點(diǎn): 其一, 由于其求解過(guò)程局部性強(qiáng), 求解過(guò)程簡(jiǎn)單 (無(wú)需求解微分方程) , 因此求解程序并行化容易, 并行效率高(Aidun & Clausen 2010). 其二,由于LBM 求解器是基于介觀(guān)尺度的方法, 因此它不僅能夠處理連續(xù)介質(zhì)問(wèn)題, 對(duì)于稀薄流動(dòng)等連續(xù)介質(zhì)假設(shè)失效的問(wèn)題也可以進(jìn)行模擬. 其三, 對(duì)于邊界條件的處理也較容易, 例如可以使用回彈方法處理無(wú)滑移壁面(Ladd 1994a, 1994b), 使用Lees-Edwards 邊界條件方法處理無(wú)界域內(nèi)的剪切流動(dòng)邊界條件(Lees & Edwards 1972)等. 此外, 對(duì)于復(fù)雜流動(dòng)問(wèn)題 (如湍流、帶界面的流動(dòng)) 都可以使用基于LBM 的求解器進(jìn)行處理(Aidun & Clausen 2010).

對(duì)于全分辨顆粒兩相流問(wèn)題, LBM 的回彈方法就可以直接地處理顆粒-流體邊界處的無(wú)滑移條件(Ladd 1994a, 1994b), 并且可以通過(guò)計(jì)算碰撞過(guò)程中的動(dòng)量交換求出顆粒受到的流體力和力矩(即FH,TH) , 從而用以推進(jìn)顆粒運(yùn)動(dòng)的求解. 但是使用這種方法處理的顆粒邊界呈臺(tái)階狀(Tenneti & Subramaniam 2014), 并存在質(zhì)量守恒無(wú)法滿(mǎn)足、動(dòng)邊界力振蕩以及伽利略不變性無(wú)法保證的問(wèn)題(Aidun & Clausen 2010). 關(guān)于這些問(wèn)題的討論以及改進(jìn)可以參看相關(guān)文獻(xiàn)(Aidun &Clausen 2010, Peng et al. 2016).

圖6格子結(jié)構(gòu)與速度的示意圖. (a) D2Q9, (b) D3Q19

2.3.3 浸沒(méi)邊界與格子玻爾茲曼混合方法

除了使用回彈模式處理顆粒-流體的邊界之外, 也可以將LBM 流體求解器和浸沒(méi)邊界法相結(jié)合 (即LBM-IBM 方法) , 通過(guò)在顆粒表面處流體施加體積力的方法來(lái)滿(mǎn)足無(wú)滑移邊界條件.這種方法的好處在于既利用了LBM 求解器的易編程、易并行化的優(yōu)點(diǎn), 同時(shí)能夠避免使用回彈方法帶來(lái)的質(zhì)量不守恒以及邊界動(dòng)量力振蕩的問(wèn)題(Aidun & Clausen 2010). 因此這種方法也被廣泛應(yīng)用于顆粒兩相流全分辨數(shù)值模擬之中(Derksen 2011, Eshghinejadfard et al. 2016).

2.3.4 顆粒碰撞

對(duì)于全分辨的顆粒兩相流模擬, 由于顆粒的尺寸不可忽略, 通常還要考慮顆粒之間的碰撞,即四向耦合. 當(dāng)不考慮流體相時(shí), 已經(jīng)有相關(guān)綜述文章詳細(xì)討論了非球形顆粒的碰撞問(wèn)題 (稱(chēng)為干碰撞) , 并且對(duì)其中涉及到的接觸檢測(cè)和碰撞模型等重點(diǎn)問(wèn)題進(jìn)行了闡述(Zhong et al. 2016).但是當(dāng)有流體存在時(shí), 除了碰撞力, 顆粒間距很小時(shí)還存在潤(rùn)滑力作用. 即便是全分辨的數(shù)值模擬也無(wú)法分辨顆粒距離很近時(shí)的潤(rùn)滑力作用, 因此有流體存在時(shí)的顆粒碰撞 (稱(chēng)為濕碰撞) 更加復(fù)雜. 一般有三種思路對(duì)濕碰撞過(guò)程進(jìn)行數(shù)值模擬. 最簡(jiǎn)單的是用簡(jiǎn)化的排斥力來(lái)模擬顆粒間的碰撞和相互作用Fip,j(Glowinski et al. 1999), 對(duì)于球形顆粒表達(dá)式為

式中?是碰撞力產(chǎn)生的閾值,?p ?1 決定了碰撞力的大小, 這兩個(gè)量均需要人為給定.di,j為兩顆粒中心位置的距離,Ri與Rj分別為兩顆粒半徑,Xi與Xj分別為兩顆粒的中心位置. 這種碰撞模型形式簡(jiǎn)單, 但是只能考慮法向碰撞力, 而且碰撞過(guò)程沒(méi)有能量損失. 第二種思路是潤(rùn)滑力模型加硬球碰撞模型(Derksen 2011, Jain et al. 2019). 其中的潤(rùn)滑力模型是基于Jeffrey 對(duì)球形顆粒間潤(rùn)滑力的理論推導(dǎo)(Jeffrey 1982), 在顆粒間距小于一定范圍 (通常是1 個(gè)歐拉網(wǎng)格) 時(shí)施加.而當(dāng)顆粒相互接觸之后, 再通過(guò)碰撞前兩顆粒的速度、碰撞過(guò)程的動(dòng)量守恒關(guān)系、和先驗(yàn)給定的恢復(fù)系數(shù)確定兩顆粒碰后的速度大小. 這種方法不需要分辨顆粒接觸后的碰撞過(guò)程, 理論上可以考慮法向和切向的動(dòng)量交換, 但是對(duì)于多體碰撞 (三個(gè)或三個(gè)以上顆粒間的碰撞) 難以處理.最后一種思路是潤(rùn)滑力加軟球碰撞模型(Xia et al. 2020). 這種思路同樣在顆粒間距很近時(shí)首先施加潤(rùn)滑力作用, 而當(dāng)顆粒產(chǎn)生相互接觸后使用基于彈性力學(xué)理論推導(dǎo)得到的碰撞力 (例如經(jīng)典的Hertz 模型) 來(lái)分辨碰撞過(guò)程. 這種方法與物理真實(shí)情況最為接近, 可以考慮碰撞過(guò)程的切向動(dòng)量交換, 也可以考慮多體碰撞問(wèn)題, 但是會(huì)由于碰撞過(guò)程剛度系數(shù)過(guò)大影響數(shù)值穩(wěn)定性, 需要很小的時(shí)間步或用多時(shí)間步的方法來(lái)處理碰撞過(guò)程. 所以有些學(xué)者通過(guò)一種“人工延長(zhǎng)碰撞時(shí)間”的方法減小碰撞過(guò)程的剛度, 從而改善碰撞模擬的數(shù)值穩(wěn)定性(Kempe & Froehlich 2012,Costa et al. 2015).

3 非球形顆粒在靜止流場(chǎng)沉降或均勻來(lái)流中轉(zhuǎn)動(dòng)

非球形顆粒在均勻來(lái)流中受到的力矩與其方位角的關(guān)系, 將決定其在均勻來(lái)流中的轉(zhuǎn)動(dòng)方式, 進(jìn)而影響其在流體中的取向和受力. 影響非球形顆粒在均勻來(lái)流下的轉(zhuǎn)動(dòng)的關(guān)鍵無(wú)量綱數(shù)是顆粒雷諾數(shù) . 本節(jié)將基于考慮不同Re的情況對(duì)該問(wèn)題進(jìn)行討論, 并且以單個(gè)顆粒在靜止流場(chǎng)中的沉降問(wèn)題為例解釋該問(wèn)題的物理意義.

Re=U∞D(zhuǎn)/ν

3.1 Stokes 流動(dòng)情況

當(dāng)忽略流體慣性效應(yīng)時(shí), 根據(jù)Jeffery 的理論(Jeffery 1922), 在均勻來(lái)流中無(wú)論顆粒的方位角如何, 都不會(huì)受到流體力矩的作用, 也就是說(shuō)顆粒會(huì)保持方位角不變. 根據(jù)這樣的結(jié)論, 對(duì)于一個(gè)非球形顆粒, 以任意的初始方位角釋放, 都會(huì)保持該角度不變. 而再根據(jù)非球形顆粒的受力公式(7), 就可以推導(dǎo)出顆粒以不同方位角沉降時(shí)的沉降速度(Ardekani et al. 2017a)

式中eg,n分別是重力方向和顆粒的回轉(zhuǎn)軸方向單位矢量, 而v1,v3分別是回轉(zhuǎn)軸與重力垂直和回轉(zhuǎn)軸與重力平行時(shí)顆粒的沉降速度, 即對(duì)于細(xì)長(zhǎng)顆粒v1對(duì)應(yīng)最大阻力沉降速度,v3對(duì)應(yīng)最小阻力沉降速度; 對(duì)于扁平顆粒則v1對(duì)應(yīng)最小阻力沉降速度,v3對(duì)應(yīng)最大阻力沉降速度. 也就是說(shuō),顆粒的沉降方位角、沉降速度矢量均可以不和重力方向重合, 顆粒以不同的取向沉降時(shí)也會(huì)有不同的沉降速度. Candelier 和Mehlig (2016)就利用細(xì)長(zhǎng)體理論得到的顆粒受力表達(dá)式, 推導(dǎo)出一根細(xì)針 (如圖7) 受重力作用沉降時(shí), 速度矢量方向β和細(xì)針的取向角α的關(guān)系

3.2 有限小Re 情況

圖7細(xì)針以任意取向在流體中沉降示意圖

在考慮流體慣性效應(yīng)后, 結(jié)果會(huì)有顯著的不同. Khyat 和Cox(1989)基于細(xì)長(zhǎng)體近似理論,利用小參數(shù)漸進(jìn)展開(kāi)方法, 推導(dǎo)出了細(xì)長(zhǎng)體 (κ=a/l ?1 ) 在均勻來(lái)流中, 當(dāng)κRe ?1 時(shí), 受到的力矩與取向角的關(guān)系. 隨后Dabade 等(2015)給出了軸對(duì)稱(chēng)橢球顆粒的流體慣性力矩在有限小Re下 (Re ?1 ) 的解析表達(dá)式式(23). 根據(jù)流體慣性力矩的形式與方位角的關(guān)系, 利用顆粒轉(zhuǎn)動(dòng)的穩(wěn)定性分析, 可以推導(dǎo)放置在均勻來(lái)流中的顆粒, 在流體慣性力矩的作用下, 最終到達(dá)的平衡狀態(tài)一定是阻力最大的取向方向.

對(duì)于沉降問(wèn)題來(lái)說(shuō), 上述結(jié)論的推論就是, 以任意初始方位角釋放的非球形顆粒, 即便在Re ?1的情況下, 也不會(huì)像完全忽略流體慣性效應(yīng)的Stokes 情況一樣保持其方位角不變沉降,而是會(huì)在流體慣性的作用下逐漸改變姿態(tài)直至穩(wěn)定到阻力最大的取向 (即細(xì)長(zhǎng)顆?；剞D(zhuǎn)軸和重力方向垂直, 扁平顆粒回轉(zhuǎn)軸和重力方向平行) , 并且沉降速度矢量也會(huì)和重力方向平行, 沉降速度大小即前面所述的最大阻力沉降速度v1.

3.3 有限大Re 情況

盡管力矩系數(shù)的擬合公式定量上誤差較大, 但是非球形顆粒在均勻來(lái)流中受到的力矩與其方位角的定性關(guān)系是明確的, 即隨著回轉(zhuǎn)軸和來(lái)流方向的夾角從0°增大到90°, 顆粒受到的力矩會(huì)先增大后減小. 并且如果顆?？梢宰杂赊D(zhuǎn)動(dòng), 最終穩(wěn)定的取向?qū)⑹亲枇ψ畲蟮姆较?Ardekani et al. 2016). 當(dāng)然如果雷諾數(shù)足夠高, 使得顆粒后方出現(xiàn)周期性的渦脫落, 顆粒方位角將在其穩(wěn)定取向方向附近振蕩.

4 非球形顆粒在線(xiàn)性剪切流中的運(yùn)動(dòng)

4.1 顆粒中心處滑移速度為零

該部分的研究最早可以追溯到1922 年, Jeffery (1922)首先在Stokes 流中研究了橢球顆粒在簡(jiǎn)單線(xiàn)性剪切作用下的轉(zhuǎn)動(dòng)行為. Jeffery 認(rèn)為忽略流體與顆粒的慣性時(shí), 顆粒是處于不受外力的狀態(tài). 此時(shí), 軸對(duì)稱(chēng)橢球顆粒在剪切流中的轉(zhuǎn)動(dòng)軌跡為穩(wěn)定的封閉曲線(xiàn), 又稱(chēng)Jeffery 軌跡. 該軌跡的形態(tài)與初始的顆粒取向有關(guān). 如果橢球顆粒的三個(gè)主軸的長(zhǎng)度均不相等, 橢球顆粒在線(xiàn)性剪切流中會(huì)呈現(xiàn)出復(fù)雜的運(yùn)動(dòng)狀態(tài), 而且在特定的形狀參數(shù)下會(huì)呈現(xiàn)混沌的轉(zhuǎn)動(dòng)的行為(Hinch &Leal, 1979, Yarin et al. 1997), 如圖8(b)～(d). Einarsson 等(2016)通過(guò)微槽道實(shí)驗(yàn)觀(guān)察到顆粒因?yàn)檩S對(duì)稱(chēng)的缺失引起的復(fù)雜的轉(zhuǎn)動(dòng)行為. 當(dāng)考慮顆粒慣性且依然忽略流體慣性作用時(shí), Lundell 和Carlsson(2010)與Challabotla 等(2015a)將Jeffery 力矩作為顆粒的合外力矩. 他們通過(guò)求解橢球顆粒的剛體轉(zhuǎn)動(dòng)的歐拉方程獲得顆粒的轉(zhuǎn)動(dòng)行為, 他們發(fā)現(xiàn)桿狀顆粒的回轉(zhuǎn)軸傾向性地垂直于渦軸并在剪切平面內(nèi)轉(zhuǎn)動(dòng), 而碟狀顆粒的回轉(zhuǎn)軸的方向傾向朝著渦軸方向且繞著回轉(zhuǎn)軸在剪切平面內(nèi)轉(zhuǎn)動(dòng), 且顆粒的轉(zhuǎn)動(dòng)形態(tài)隨顆粒慣性的不同而呈現(xiàn)不同的取向軌跡線(xiàn) (Lundell &Carlsson 2010, Challabotla et al. 2015a) 如圖8(e)和圖8(f). 對(duì)于三軸不等顆粒, Lundell(2011b)發(fā)現(xiàn)顆粒慣性會(huì)讓顆粒從混沌的轉(zhuǎn)動(dòng)的狀態(tài)轉(zhuǎn)變?yōu)轭w粒繞自身最短軸在剪切平面內(nèi)周期轉(zhuǎn)動(dòng)的狀態(tài), 如圖8(g). 為了進(jìn)一步解釋顆粒慣性與形狀對(duì)顆粒行為的作用, Yarin 等(1997)最早使用了線(xiàn)性穩(wěn)定性分析中的Floquet 理論研究了顆粒形狀對(duì)轉(zhuǎn)動(dòng)行為的影響, 并發(fā)現(xiàn)趨近桿狀的三軸不等顆粒的運(yùn)動(dòng)行為更容易出現(xiàn)不穩(wěn)定轉(zhuǎn)動(dòng)的行為. Lundell (2011b)也通過(guò)Floquet 分析理論給出顆粒慣性對(duì)三軸不等顆粒的轉(zhuǎn)動(dòng)穩(wěn)定性的影響. 而Einarsson 等(2014)對(duì)顆粒方程基于弱顆粒慣性進(jìn)行攝動(dòng)展開(kāi), 并利用Floquet 理論解析地得到軸對(duì)稱(chēng)橢球顆粒在弱顆粒慣性下自旋與翻轉(zhuǎn)模態(tài)運(yùn)動(dòng)的穩(wěn)定性. Cui Z 等(2019)系統(tǒng)地研究了軸對(duì)稱(chēng)與三軸不等橢球顆粒繞各個(gè)顆粒主軸方向轉(zhuǎn)動(dòng)的穩(wěn)定性, 并從理論上提供了更加廣泛的顆粒慣性與形狀的參數(shù)下顆粒不同轉(zhuǎn)動(dòng)模態(tài)的穩(wěn)定性. 與此同時(shí), 當(dāng)顆粒在簡(jiǎn)單線(xiàn)性剪切平面中進(jìn)行準(zhǔn)二維轉(zhuǎn)動(dòng)時(shí), 不管顆粒形狀與慣性如何變化, 顆粒均是周期轉(zhuǎn)動(dòng). Lundell 和Carlsson (2010)發(fā)現(xiàn)顆粒的轉(zhuǎn)動(dòng)周期在顆粒慣性的臨界值附近存在較大的階躍, 當(dāng)顆粒慣性遠(yuǎn)小于臨界值時(shí), 顆粒的轉(zhuǎn)動(dòng)周期只與形狀有關(guān)但與慣性無(wú)關(guān), 而在大于臨界值時(shí), 不同形狀的顆粒轉(zhuǎn)動(dòng)周期均趨于相同值. 當(dāng)流體的剪切形式隨時(shí)間進(jìn)行周期變化, 在特定的顆粒形狀與慣性范圍內(nèi), 軸對(duì)稱(chēng)橢球顆粒的轉(zhuǎn)動(dòng)也會(huì)呈現(xiàn)出混沌的運(yùn)動(dòng)現(xiàn)象(Nilsen et al. 2013).

圖8(a)顆粒在簡(jiǎn)單剪切流示意圖; (b)～(g)橢球顆粒在剪切流中的轉(zhuǎn)動(dòng)時(shí)代表橢球顆粒方向的軌跡圖. 其中(b)～(d)無(wú)慣性橢球顆粒, (e)～(g)以最長(zhǎng)軸定義的顆粒轉(zhuǎn)動(dòng)慣性分布為 S tr =669, 646, 100 ;形 狀 參 數(shù)(b)(e) λ =5 , (c)(f) λ =1/5 , (d)(g)ka =1,kb =10-0.65,kc =0.1

如果考慮流體慣性的影響, 非球形顆粒的轉(zhuǎn)動(dòng)將會(huì)變得更加復(fù)雜. Qi 和Luo (2002, 2003)較早利用格子玻爾茲曼方法研究了有限尺寸橢球顆粒在剪切流中的轉(zhuǎn)動(dòng)行為, 并發(fā)現(xiàn)隨著顆粒雷諾數(shù)的增大以及顆粒形狀的變化, 顆粒呈現(xiàn)出多種轉(zhuǎn)動(dòng)與取向模態(tài). 隨后, Yu 等(2007a), Huang等(2012) 以及Rosén 等 (2014) 進(jìn)一步研究了軸對(duì)稱(chēng)橢球顆粒在剪切流中的行為, 發(fā)現(xiàn)并總結(jié)出在不同顆粒剪切雷諾數(shù)下橢球顆粒的不同運(yùn)動(dòng)模態(tài) (如圖9) , 隨著顆粒慣性的增強(qiáng), 桿狀與碟狀顆粒的運(yùn)動(dòng)模態(tài)發(fā)生變化. 表2 總結(jié)了長(zhǎng)寬比為λ=2 和4 時(shí)桿狀顆粒隨顆粒剪切雷諾數(shù)變化時(shí)不同的運(yùn)動(dòng)模態(tài). 表3 總結(jié)了長(zhǎng)寬比λ=1/2 的碟狀顆粒在不同雷諾數(shù)下的轉(zhuǎn)動(dòng)模態(tài). Rosén 與其合作者(2015a, 2015b, 2016, 2017a, 2017b)進(jìn)一步從定量角度深入探討有限尺寸橢球顆粒在剪切流中的復(fù)雜行為 (包括不同模態(tài)以及模態(tài)之間的轉(zhuǎn)換條件) , 并給出了定量研究的網(wǎng)格要求.他們發(fā)現(xiàn)顆粒慣性與流體慣性對(duì)顆粒轉(zhuǎn)動(dòng)行為的影響呈現(xiàn)出競(jìng)爭(zhēng)的關(guān)系, 進(jìn)而導(dǎo)致相同形狀顆粒在不同初始位置下可能存在不同的轉(zhuǎn)動(dòng)行為. 同時(shí), Rosén (2017a)發(fā)現(xiàn)了具有較大長(zhǎng)寬比的桿狀顆粒在流體慣性與顆粒慣性的相互作用下, 其自旋模態(tài)(回轉(zhuǎn)軸與渦軸方向平行)可能會(huì)出現(xiàn)Shilnikov 分岔甚至混沌的行為. 當(dāng)顆粒為三軸不等橢球顆粒時(shí), 在極低顆粒雷諾數(shù)下三軸不等顆?？赡軙?huì)出現(xiàn)混沌行為(Hinch & Leal 1979, Yarin et al. 1997, Lundell 2011, Cui Z et al.2019), Rosén 等 (2017b)較系統(tǒng)地研究了流體慣性對(duì)三軸不等橢球顆粒在剪切流中的運(yùn)動(dòng)行為.研究表明在流體慣性作用下三軸不等橢球顆粒的行為與極低雷諾數(shù)的情況下類(lèi)似, 不過(guò)流體慣性會(huì)對(duì)顆粒繞中間軸的轉(zhuǎn)動(dòng)有一定穩(wěn)定作用. 在理論研究方面, Meibohm 等(2016)通過(guò)攝動(dòng)理論得到了軸對(duì)稱(chēng)橢球顆粒在自旋時(shí)的角速度隨顆粒剪切雷諾數(shù)的關(guān)系. 而Dabade 等(2016)與Einarsson 等(2015)利用互等原理系統(tǒng)地研究在小剪切雷諾數(shù)下流體慣性對(duì)顆粒轉(zhuǎn)動(dòng)行為的影響.

圖9顆粒在剪切流中轉(zhuǎn)動(dòng)的不同模態(tài)(Rosén et al. 2014)

表2 有限尺寸桿狀顆粒在剪切流中不同 R es 對(duì)應(yīng)的轉(zhuǎn)動(dòng)模態(tài)

表3 有限尺寸碟狀顆粒在剪切流中不同 R es 對(duì)應(yīng)的轉(zhuǎn)動(dòng)模態(tài)

在Jeffery (1922)推導(dǎo)出橢球顆粒在Stokes 流中的力矩方程的同時(shí), 也提出一個(gè)關(guān)于顆粒最終穩(wěn)定轉(zhuǎn)動(dòng)狀態(tài)應(yīng)該處于系統(tǒng)最小能量耗散的猜想. Jeffery 認(rèn)為桿狀顆粒的回轉(zhuǎn)軸平行于渦軸進(jìn)行自旋和碟狀顆粒的回轉(zhuǎn)軸垂直于渦軸在速度梯度平面翻轉(zhuǎn)為最終的穩(wěn)定狀態(tài). 盡管Taylor(1923)通過(guò)實(shí)驗(yàn)觀(guān)察顆粒轉(zhuǎn)動(dòng)的軌跡驗(yàn)證了Jeffery 提出的最小能量耗散的猜想, 但在近些年的直接數(shù)值模擬的結(jié)果中均不支持該猜想(Qi & Luo 2003, Rosén et al. 2014, Huang et al. 2012,Mao & Alexeev 2014). Huang 等 (2012) 的直接數(shù)值模擬的結(jié)果表明, 桿狀顆粒的回轉(zhuǎn)軸在速度梯度平面的翻轉(zhuǎn)和碟狀顆粒的回轉(zhuǎn)軸平行于渦軸進(jìn)行自旋是最終穩(wěn)定的狀態(tài), 他們對(duì)應(yīng)的是最大能量耗散的狀態(tài). Mao 和Alexeev (2014) 則發(fā)現(xiàn)對(duì)于長(zhǎng)寬比為2 的桿狀顆粒, 當(dāng)顆粒剪切雷諾數(shù)大于60, 顆粒的轉(zhuǎn)動(dòng)會(huì)從翻轉(zhuǎn)變?yōu)樽孕B(tài). Einarsson 等 (2015) 通過(guò)攝動(dòng)理論對(duì)顆粒運(yùn)動(dòng)方程基于小剪切雷諾數(shù)進(jìn)行展開(kāi), 理論表明非常扁平的碟狀顆粒的翻轉(zhuǎn)與自旋在小剪切雷諾數(shù)下都可能是穩(wěn)定的模態(tài). 這些研究成果意味著Jeffery 提出的最小能量耗散的猜想可能存在一定的適用范圍, 比如顆粒形狀范圍與顆粒剪切雷諾數(shù)范圍等.

4.2 顆粒中心處滑移速度非零

對(duì)于大多數(shù)實(shí)際情況, 顆粒在流場(chǎng)中輸運(yùn)過(guò)程中一定會(huì)存在平動(dòng)的滑移速度. 僅具有滑移速度的顆粒問(wèn)題已經(jīng)被大量的研究, 相關(guān)內(nèi)容可以參見(jiàn)第3 節(jié). 但同時(shí)考慮剪切與滑移的研究目前還缺少系統(tǒng)研究. Li 等(2019)研究了將顆粒中心固定在僅有一個(gè)壁面滑動(dòng)的Couette 槽道流動(dòng)的中心位置的轉(zhuǎn)動(dòng)情況. 其結(jié)果展示了橢球顆粒在平動(dòng)滑移產(chǎn)生的流體慣性力矩與剪切產(chǎn)生的力矩作用下產(chǎn)生的轉(zhuǎn)動(dòng)行為, 但尚未對(duì)其運(yùn)動(dòng)行為展開(kāi)系統(tǒng)與細(xì)致的深入研究工作(圖10).

5 非球形顆粒在均勻各向同性湍流中的運(yùn)動(dòng)

均勻各向同性湍流是一種最簡(jiǎn)單的湍流場(chǎng), 其各階湍流統(tǒng)計(jì)特性在空間處處相等且不隨坐標(biāo)系的剛性轉(zhuǎn)動(dòng)而改變(張兆順等 2017). 為了研究非球形顆粒與湍流的相互作用, 均勻各向同性湍流是一個(gè)比較合適的選擇. 近些年, 有關(guān)非球形顆粒在均勻各向同性湍流中運(yùn)動(dòng)的研究比較活躍, 接下來(lái)將主要針對(duì)點(diǎn)顆粒與全分辨顆粒方法, 對(duì)目前的研究現(xiàn)狀進(jìn)行總結(jié).

5.1 點(diǎn)顆粒橢球顆粒

5.1.1 無(wú)慣性橢球顆粒

圖10顆粒中心處滑移速度非零時(shí)的剪切流模型示意圖

對(duì)于顆粒尺寸遠(yuǎn)小于湍流的Kolmogorov 尺寸, 可忽略流體慣性以及顆粒慣性的作用. 此時(shí)顆粒軌跡可近似跟隨流體質(zhì)點(diǎn), 顆粒的轉(zhuǎn)動(dòng)可以根據(jù)當(dāng)?shù)氐牧黧w速度梯度和顆粒取向信息通過(guò)Jeffery 方程求得, 而顆粒取向需要沿顆粒軌跡進(jìn)行積分求得. 無(wú)慣性橢球顆粒因?yàn)槠淠Ｐ秃?jiǎn)單并且具有很好的跟隨性而被廣泛應(yīng)用在非球形顆粒湍流的理論研究中. 盡管實(shí)際問(wèn)題中顆粒慣性和流體慣性并不能完美與無(wú)慣性顆粒的假設(shè)符合, 但大量的實(shí)驗(yàn)觀(guān)測(cè)結(jié)果表明, 無(wú)慣性橢球顆粒模型能較好預(yù)測(cè)微小顆粒(流體慣性和顆粒慣性影響非常小)的復(fù)雜行為 (Parsa et al. 2012,Marcus et al. 2014, Fries et al. 2017, Einarsson et al. 2016). Pumir 和Wilkinson (2011)和Ni 等(2014)報(bào)道了桿狀顆粒的取向與流體的渦量方向具有強(qiáng)相關(guān). 同時(shí), 他們還發(fā)現(xiàn)桿狀顆粒沿變形率張量的最大拉伸方向相關(guān)較弱, 但與變形率張量的第二特征方向相關(guān)反而較強(qiáng). 但是顆粒取向與渦量的相關(guān)性主要在顆粒尺寸小于耗散尺度時(shí)較強(qiáng). 如果顆粒尺寸超過(guò)流體的耗散尺度, 且在流體的慣性尺度的范疇, 那么桿狀顆粒的取向是與變形率張量的最大拉伸方向相關(guān)最強(qiáng)(Pujara et al. 2019). 而對(duì)于碟狀顆粒, 其回轉(zhuǎn)軸傾向垂直于渦量方向, 且與變形率的最大壓縮方向具有較強(qiáng)相關(guān) (Gustavsson et al. 2014, Byron et al. 2015) . Ni 等(2014)發(fā)現(xiàn)桿狀顆粒的取向與流體拉格朗日拉伸方向具有非常強(qiáng)的相關(guān)性. 他們認(rèn)為長(zhǎng)寬比大于10 的桿狀顆粒與流體拉格朗日拉伸方向具有非常好的一致性. 與此同時(shí), Zhao 等(2019b)發(fā)現(xiàn)顆粒的取向在空間的分布并不連續(xù), 這意味著兩個(gè)距離很近(小于Kolmogorov 尺寸)的顆粒依然存在大的夾角. 在顆粒轉(zhuǎn)動(dòng)方面, Parsa 等(2012)、Marcus 等(2014)和Byron 等(2015)通過(guò)實(shí)驗(yàn)和數(shù)值模擬得到了均勻各向同性湍流中顆粒翻轉(zhuǎn)率(tumbling)與自旋率(spinning)隨顆粒形狀變化的函數(shù)關(guān)系 (如圖11(a)所示) . 該圖表明碟狀顆粒的翻轉(zhuǎn)率強(qiáng)于桿狀顆粒, 但碟狀顆粒的自旋率弱于桿狀顆粒. 而當(dāng)顆粒長(zhǎng)寬比大于3 和小于1/3 時(shí)顆粒轉(zhuǎn)動(dòng)變化不明顯. 然而, 在關(guān)于顆粒轉(zhuǎn)動(dòng)的前期的統(tǒng)計(jì)理論模型的研究中, 顆粒的轉(zhuǎn)動(dòng)行為對(duì)于碟狀與桿狀顆粒沒(méi)有區(qū)別, 其轉(zhuǎn)動(dòng)率隨顆粒形狀的變化應(yīng)該是基于λ=1 對(duì)稱(chēng)的(Shin & Koch 2005, Byron et al. 2015). 在該模型中, 湍流脈動(dòng)隨機(jī)且顆粒取向隨機(jī), 顆粒取向與流動(dòng)相互獨(dú)立. 而Parsa 等(2012)的工作從實(shí)驗(yàn)與數(shù)值上均說(shuō)明了顆粒取向與流動(dòng)結(jié)構(gòu)的相關(guān)的重要性. Gustavsson 等 (2014) 從理論上解釋在傳統(tǒng)基于隨機(jī)湍流的模型預(yù)測(cè)的顆粒行為與真實(shí)湍流中顆粒轉(zhuǎn)動(dòng)行為差異的原因. 他們認(rèn)為碟狀顆粒的持續(xù)翻轉(zhuǎn)與顆粒軌跡經(jīng)過(guò)高渦量區(qū)有關(guān), 而傳統(tǒng)的模型還不能捕捉. Pujara 等(2021)研究不同形狀的顆粒的轉(zhuǎn)動(dòng)隨著對(duì)直接數(shù)值模擬的湍流場(chǎng)過(guò)濾尺寸變化的關(guān)系. 研究結(jié)果表明無(wú)論用Kolmogorov 尺寸還是用積分尺度過(guò)濾流場(chǎng), 顆粒的轉(zhuǎn)動(dòng)行為與顆粒形狀的關(guān)系與過(guò)濾尺寸影響不大, 這表明非球形顆粒的轉(zhuǎn)動(dòng)可能體現(xiàn)了一種拉格朗日尺度的不變性. Chen 等(2016)發(fā)現(xiàn)目前的大渦模擬方法中流體相的亞格子模式對(duì)非球形顆粒的轉(zhuǎn)動(dòng)行為的影響較大, 使得其結(jié)構(gòu)與直接數(shù)值模擬的結(jié)果具有較大偏差. 近日, Pujara 等(2021)提出的不同尺寸過(guò)濾后的拉格朗日尺度不變性可能會(huì)對(duì)非球形顆粒的亞格子模型的建立具有參考意義. 如果顆粒尺寸大于Kolmogorov 尺寸, Pujara 等(2019)發(fā)現(xiàn)桿狀顆粒的翻轉(zhuǎn)率與顆粒自身的長(zhǎng)度滿(mǎn)足-4/3 的冪律關(guān)系. 對(duì)于三軸不等顆粒, Chevillard 和Meneveau (2013)研究了三軸不等顆粒在各向同性湍流中的取向行為并與隨機(jī)模型進(jìn)行對(duì)比. Pujara & Variano (2017)研究了三軸不等橢球顆粒在各向同性湍流中的轉(zhuǎn)動(dòng)行為. 研究結(jié)果表明顆粒的擬渦能主要被分配給最長(zhǎng)軸, 并且顆粒繞著其最長(zhǎng)軸的轉(zhuǎn)動(dòng)是最持久的. 同時(shí), 研究發(fā)現(xiàn)湍流中渦量-變形率的相關(guān)使得顆粒的總擬渦能與顆粒形狀的關(guān)系較弱.

圖11無(wú)慣性橢球顆粒的(a)轉(zhuǎn)動(dòng)率隨形狀變化的關(guān)系(數(shù)值模擬數(shù)據(jù)來(lái)自Byron 等 (2015), 兩處實(shí)驗(yàn)數(shù)據(jù)分別來(lái)自Marcus 等 (2014)與Parsa 等(2012))與(b)在湍流中運(yùn)動(dòng)示意圖

5.1.2 慣性橢球顆粒

當(dāng)考慮顆粒慣性的作用時(shí), 近期的實(shí)驗(yàn)發(fā)現(xiàn), 纖維自身的慣性會(huì)降低纖維自身的轉(zhuǎn)動(dòng)(Bounoua et al. 2018). 然而, 非球形慣性顆粒在各向同性湍流的數(shù)值模擬研究卻比較少. Roy 等(2018)僅研究了顆粒平動(dòng)慣性對(duì)顆粒運(yùn)動(dòng)的影響, 卻忽略了顆粒轉(zhuǎn)動(dòng)慣性的作用. 他們認(rèn)為在各向同性湍流中非球形顆粒的轉(zhuǎn)動(dòng)弛豫時(shí)間要明顯小于平動(dòng)的弛豫時(shí)間. Zhao 等(2015)在槽道湍流的中部(流動(dòng)趨近于各向同性湍流)研究了慣性非球形顆粒的轉(zhuǎn)動(dòng)行為隨顆粒慣性的變化. 其結(jié)果也展示了顆粒慣性會(huì)削弱顆粒自身的轉(zhuǎn)動(dòng)行為. Gustavsson 等(2014)則是通過(guò)建立隨機(jī)湍流模型研究了慣性非球形顆粒的轉(zhuǎn)動(dòng)行為. 關(guān)于顆粒慣性對(duì)顆粒在各向同性湍流中的取向與轉(zhuǎn)動(dòng)行為的影響依然缺少系統(tǒng)與深入的研究.

5.1.3 橢球顆粒沉降

在早期關(guān)于球形顆粒的沉降問(wèn)題中, 湍流會(huì)加速顆粒的沉降(Maxey 1987, Wang & Maxey 1993, Bec et al. 2014). 這是因?yàn)橹仡w粒會(huì)傾向性地聚集在流動(dòng)向下的區(qū)域(Wang & Maxey 1993, Bec et al. 2014). 對(duì)于非球形顆粒而言, 顆粒沉降的過(guò)程比球形顆粒的沉降更加復(fù)雜. 一方面是形狀的各向異性導(dǎo)致沉降速度不僅與湍流強(qiáng)度有關(guān)還與顆粒取向有關(guān)(Siewert et al.2014a); 另一方面, 非球形顆粒的碰撞頻率要比球形顆粒更加強(qiáng)烈(Siewert et al. 2014b). 為忽略顆粒慣性的影響, Ardekani 等(2017a)直接在無(wú)慣性顆粒的模型上添加了重力沉降速度. 其研究結(jié)果表明顆粒形狀的增加會(huì)影響顆粒聚集以及平均的沉降速度. 其背后可能的原因有兩個(gè): 一個(gè)是形狀對(duì)顆粒轉(zhuǎn)動(dòng)的影響, 另一個(gè)是提高了阻力的各向異性, 而前者的作用更加明顯. Jucha 等(2018) 通過(guò)研究碟狀顆粒的沉降過(guò)程時(shí)發(fā)現(xiàn)在湍流能量耗散率較小時(shí)顆粒的碰撞主要由相鄰顆粒的取向差異引起的重力沉降導(dǎo)致. 隨著湍流能量耗散率的增加, 湍流的脈動(dòng)主要誘導(dǎo)了顆粒的碰撞行為. Gustavsson 等(2017)利用簡(jiǎn)單的高斯隨機(jī)模型就準(zhǔn)確預(yù)測(cè)了非球形重顆粒在湍流中的取向行為.

隨著Dabade 等(2015)推導(dǎo)出軸對(duì)稱(chēng)橢球顆粒流體慣性力矩的解析表達(dá)式, Gustavsson 等(2019)、Sheikh 等(2020) 和Anand 等(2020)意識(shí)到流體慣性力矩在非球形顆粒在湍流中沉降的重要性. 流體慣性力矩主要由顆粒與流場(chǎng)的滑移速度誘導(dǎo), Gustavsson 等(2019)發(fā)現(xiàn)當(dāng)Rep ?1時(shí), 流體慣性力矩的作用依然不可忽略. Sheikh 等(2020)引入無(wú)量綱參數(shù) ? 衡量流體慣性力矩與黏性力矩. 當(dāng) ??1 時(shí), 顆粒傾向朝著阻力最大的方向進(jìn)行沉降.

5.2 有限尺寸顆粒

而在各向同性湍流中進(jìn)行全分辨顆粒的模擬的研究到目前為止依然比較少. Schneiders 等(2017, 2019)在各向湍流中研究了尺寸接近流體Kolmogorov 尺寸的非球形顆粒運(yùn)動(dòng)行為以及顆粒與流體之間的相互作用, 并與傳統(tǒng)的拉格朗日雙向耦合點(diǎn)顆粒模型進(jìn)行比較(圖12). 結(jié)果表明非球形點(diǎn)顆粒模型在預(yù)測(cè)取向行為上的準(zhǔn)確性, 但在顆粒速度、顆粒與流場(chǎng)的相互作用等方面依然還存在明顯差距. 如何基于全分辨顆粒的研究結(jié)果去改進(jìn)拉格朗日點(diǎn)顆粒模型在未來(lái)依然是一個(gè)值得研究的方向.

6 非球形顆粒在壁湍流中的運(yùn)動(dòng)

在實(shí)際工程問(wèn)題中, 湍流具有各向異性且在空間上分布不均勻. 壁湍流就是一種典型的強(qiáng)剪切和強(qiáng)各向異性的湍流模型. 由于壁面的存在, 在壁面附近不僅存在強(qiáng)的速度剪切, 而且還存在非常豐富的湍流結(jié)構(gòu), 比如流向渦、條帶、猝發(fā)等(許春曉 2015). 隨著雷諾數(shù)的升高, 壁湍流的外區(qū)存在大尺度結(jié)構(gòu), 這些大尺度結(jié)構(gòu)的流向長(zhǎng)度通常為2～3 倍、13～15 倍甚至20 倍外區(qū)尺度(許春曉 2015). Zhang 等(2001)最早在壁湍流中進(jìn)行非球形顆粒兩相流的直接數(shù)值模擬研究, 隨后研究學(xué)者們陸續(xù)在壁湍流中展開(kāi)了全方面的非球形顆粒兩相流研究. 目前對(duì)非球形顆粒在壁湍流的運(yùn)動(dòng)行為特征有了基本的認(rèn)識(shí). 本節(jié)圍繞近些年關(guān)于非球形顆粒在壁湍流運(yùn)動(dòng)的研究進(jìn)展展開(kāi)討論.

6.1 點(diǎn)顆粒橢球顆粒

6.1.1 無(wú)慣性橢球顆粒

圖12有限尺寸桿狀顆粒在湍流中分布示意圖(a)與(b)自適應(yīng)網(wǎng)格加密(Schneiders et al. 2019)

由于壁面的存在, 壁面附近的流動(dòng)呈現(xiàn)出明顯的各向異性與不均勻性的特點(diǎn). Challabotla等(2015b) 最早研究了無(wú)慣性橢球顆粒在槽道湍流中的取向行為. 在近壁面附近, 桿狀顆粒的回轉(zhuǎn)軸方向傾向朝著流向方向, 而碟狀顆粒的回轉(zhuǎn)軸傾向指向壁面的法向方向(如圖13(a)(b)).但在遠(yuǎn)離壁面的槽道中部附近, 顆粒的取向行為趨于同顆粒在各向同性湍流中的行為類(lèi)似, 這是因?yàn)樵诓鄣乐胁康耐牧骶哂芯植扛飨蛲缘奶匦?Andersson et al. 2015). 同時(shí), 近壁面附近存在復(fù)雜的湍流相干結(jié)構(gòu), Cui Z 等(2021)首次觀(guān)察到非球形顆粒取向與流向渦結(jié)構(gòu)、上拋和下掃事件以及速度條帶之間的相關(guān), 并提出顆粒的取向行為在近壁面附近可以被定性地劃分為剪切主導(dǎo)區(qū)、結(jié)構(gòu)主導(dǎo)區(qū)以及各向同性區(qū)(如圖14). 在剪切主導(dǎo)區(qū), 顆粒在強(qiáng)剪切的作用下其旋轉(zhuǎn)會(huì)受到抑制, 桿狀顆粒與碟狀顆粒會(huì)在特定的朝向停留較長(zhǎng)時(shí)間 (Challabotla et al. 2015b, Zhao et al. 2015) . 逐漸遠(yuǎn)離壁面, 流場(chǎng)的剪切減弱, 湍流結(jié)構(gòu)增加, 進(jìn)入結(jié)構(gòu)主導(dǎo)區(qū)域. 此時(shí), 桿狀顆粒與碟狀顆粒均會(huì)受到結(jié)構(gòu)的影響, 進(jìn)而使得顆粒在流向渦周?chē)?、上拋與下掃事件中以及高低速條帶中的取向行為產(chǎn)生明顯差異(Cui Z et al. 2021). 該類(lèi)似現(xiàn)象在早期的實(shí)驗(yàn)中也得到了觀(guān)察, 即桿狀顆粒在高速條帶中取向相對(duì)隨機(jī)而在低速條帶取向更趨近于流向方向(Abbasi Hoseini et al. 2015). 但是該工作將此現(xiàn)象簡(jiǎn)單地歸結(jié)為顆粒與壁面碰撞的有限尺寸的影響, 并未深入分析和探討背后的原因. 而Cui Z 等(2021)的工作已經(jīng)消除了潛在的顆粒與壁面之間由于顆粒有限尺寸導(dǎo)致的碰撞效應(yīng)的影響. 當(dāng)顆粒位置進(jìn)一步遠(yuǎn)離壁面, 剪切趨近于零, 而湍流也趨于各向同性, 此時(shí)顆粒的轉(zhuǎn)動(dòng)與取向行為逐漸與其在各向同性湍流中的行為類(lèi)似(Zhao et al.2015). 不過(guò), 在中高等雷諾數(shù)槽道湍流中, 由于流場(chǎng)中存在大尺度結(jié)構(gòu)(比如等動(dòng)量區(qū)與非等動(dòng)量區(qū)), Jie 等(2019)發(fā)現(xiàn)等動(dòng)量區(qū)的內(nèi)外顆粒的轉(zhuǎn)動(dòng)統(tǒng)計(jì)行為存在明顯差異, 圖15 展示了Reτ=1000槽道湍流中桿狀顆粒的取向分布瞬時(shí)圖.

圖13橢球顆粒在槽道湍流中的取向分布. (a)(b)無(wú)慣性橢球顆粒(Challabotla et al. 2015b); (c)(d)St =30 橢球顆粒

同時(shí), 由于壁面剪切的存在, 桿狀顆粒與碟狀顆粒在近壁面附近的取向均與流體的渦量方向接近垂直 (Zhao et al. 2015) . Zhao 和Andersson (2016)發(fā)現(xiàn)桿狀顆粒與碟狀顆粒分別與流體拉格朗日的拉伸與壓縮方向有關(guān), 但相關(guān)程度并沒(méi)有其在各向同性湍流中強(qiáng). Cui Z 等(2020)進(jìn)一步研究了桿狀顆粒的取向與流體拉格朗日拉伸方向的相關(guān), 并發(fā)現(xiàn)Ni 等(2014)提出的長(zhǎng)寬比大于10 的桿狀顆粒與流體拉格朗日拉伸方向一致的結(jié)論在槽道中并不是處處成立. Cui Z 等(2020)發(fā)現(xiàn)長(zhǎng)寬比大于10 的桿狀顆粒行為與流體拉格朗日拉伸方向在黏性底層存在明顯差異.以圖16 為例, 其展示了長(zhǎng)寬比為23.3 的桿狀顆粒與流體拉格朗日拉伸方向沿同一條顆粒軌跡線(xiàn)的取向行為. 圖16 的結(jié)果顯示流體拉格朗日拉伸方向在靠近上壁面時(shí)主要朝著流向方向而桿狀顆粒在壁面附近持續(xù)翻轉(zhuǎn). Cui Z 等(2020)認(rèn)為這一現(xiàn)象主要與近壁面的速度剪切和速度梯度脈動(dòng)有關(guān). 當(dāng)考慮三軸不等顆粒時(shí), Challabotla 等(2016a)發(fā)現(xiàn)顆粒同時(shí)具有碟狀與桿狀顆粒的特性, 即顆粒像桿狀顆粒一樣繞著最長(zhǎng)軸轉(zhuǎn)動(dòng)也像碟狀顆粒一樣繞著最短軸轉(zhuǎn)動(dòng).

如果壁面剪切趨于零, Yang 等(2018)通過(guò)構(gòu)造Couette-Poiseuille 流動(dòng)發(fā)現(xiàn)無(wú)慣性的桿狀顆粒在零剪切壁面附近取向趨向隨機(jī)而碟狀顆粒依然指向壁面法向. 同時(shí), 他們的研究結(jié)果還表明顆粒的各向異性的轉(zhuǎn)動(dòng)模態(tài)主要與顆粒取向的各向異性有關(guān), 而壁面平均剪切對(duì)轉(zhuǎn)動(dòng)模態(tài)的影響起次要作用.

6.1.2 慣性橢球顆粒

圖14顆粒在近壁流向渦結(jié)構(gòu)影響下的取向行為及區(qū)域劃分(Cui Z et al. 2021). (a)(b)分別為細(xì)長(zhǎng)桿狀與扁平顆粒在瞬時(shí)流向渦附近的取向分布; (c)(d)條件系綜平均后的顆粒取向分布, 其中(c)細(xì)長(zhǎng)桿狀顆粒與流向夾角余弦值 | cosθx| ; (d) 扁平顆粒與展向夾角余弦值 | cosθz| ; (e)依據(jù)顆粒取向行為特點(diǎn)進(jìn)行的區(qū)域劃分示意圖

圖15桿狀顆粒在雷諾數(shù) R eτ =1000 槽道湍流分布圖. 其中云圖為流向速度, 白色等值線(xiàn)代表0.95 倍槽流中心平均速度

圖16無(wú)慣性桿狀顆粒在壁面取向行為與流體拉格朗日拉伸方向的差異(Cui Z et al. 2020)

在壁湍流中, 考慮顆粒慣性的研究要早于無(wú)慣性顆粒的研究, Zhang 等(2001)最早實(shí)現(xiàn)了壁湍流中慣性桿狀顆粒的輸運(yùn)與沉積特性. 隨后, Mortensen 等(2008a, 2008b)、Marchioli 等(2010,2013, 2016)、Zhao 等(2013, 2014, 2015)、Challabotla 等(2015c)和Yuan 等(2018)的研究表明慣性橢球顆粒傾向性地聚集在低速度條帶, 且與顆粒形狀的影響較小, 這與慣性球形顆粒的輸運(yùn)特性類(lèi)似. 然而顆粒的轉(zhuǎn)動(dòng)與取向行為依賴(lài)顆粒形狀. 對(duì)于桿狀顆粒, 隨著顆粒長(zhǎng)寬比的增大,桿狀顆粒在近壁面傾向朝著流向方向, 但隨著顆粒慣性的增大, 傾向性取向行為減弱, 進(jìn)而桿狀顆粒傾向在強(qiáng)剪切平面內(nèi)翻滾地向下游輸運(yùn). 對(duì)于碟狀顆粒, 隨著顆粒非球形度的增大, 弱慣性的碟狀顆粒在近壁面傾向朝著壁面法向, 但隨著顆粒慣性的增大, 碟狀顆粒像車(chē)輪一樣在壁面滾動(dòng)著向下游輸運(yùn). Zhao 等(2014)與Marchioli 等(2016)分別研究了桿狀顆粒與流場(chǎng)之間平動(dòng)滑移速度與轉(zhuǎn)動(dòng)滑移角速度隨顆粒慣性與形狀的影響. 其結(jié)果表明隨著慣性逐漸趨近于零, 顆粒與流場(chǎng)之間的轉(zhuǎn)動(dòng)滑移角速度并不像平動(dòng)滑移速度那樣趨近于零, 依然存在明顯的滑移角速度.Zhao 等(2015)對(duì)槽道湍流中慣性橢球顆粒的轉(zhuǎn)動(dòng)進(jìn)行分析, 發(fā)現(xiàn)壁湍流中的近壁轉(zhuǎn)動(dòng)模態(tài)與遠(yuǎn)離壁面的轉(zhuǎn)動(dòng)模態(tài)是完全不同的. 在遠(yuǎn)離壁面的轉(zhuǎn)動(dòng)模態(tài), 碟狀顆粒傾向繞非回轉(zhuǎn)軸垂直于渦量方向進(jìn)行翻轉(zhuǎn), 而桿狀顆粒繞回轉(zhuǎn)軸沿渦量方向自旋. 隨著顆粒慣性的增大, 這種趨勢(shì)被削弱.然而在近壁區(qū)域, 顆粒的慣性反而增加了碟狀顆粒繞回轉(zhuǎn)軸在剪切平面轉(zhuǎn)動(dòng)而桿狀顆粒繞非回轉(zhuǎn)軸轉(zhuǎn)動(dòng). 進(jìn)一步, Zhao 等(2019a)系統(tǒng)研究軸對(duì)稱(chēng)橢球顆粒在槽道不同高度時(shí)的轉(zhuǎn)動(dòng)模態(tài)隨剪切、湍流強(qiáng)度與顆粒慣性的影響, 并研究了從壁面到槽道中部?jī)煞N轉(zhuǎn)動(dòng)模態(tài)的轉(zhuǎn)變. 最近,Michel 和Arcen (2021a)的研究表明非球形顆粒的統(tǒng)計(jì)量需要較長(zhǎng)的時(shí)間才能統(tǒng)計(jì)穩(wěn)定, 同時(shí),Marchioli 等(2016)的報(bào)道表明顆粒的轉(zhuǎn)動(dòng)和取向行為比平動(dòng)的統(tǒng)計(jì)更快統(tǒng)計(jì)穩(wěn)定. Michel 和Arcen (2021b)研究了不同槽道湍流的雷諾數(shù)(從Reτ=180 到550)對(duì)慣性桿狀顆粒的行為進(jìn)行研究. 其結(jié)果表明隨著雷諾數(shù)的增大, 顆粒在壁面的傾向性取向與轉(zhuǎn)動(dòng)行為有輕微的減弱.

6.1.3 重力沉降考慮重力沉降時(shí), 目前大部分相關(guān)研究中主要考慮重力方向與槽道流向平行的情況, 即豎直槽道. 當(dāng)槽道流動(dòng)方向與重力方向一致, 稱(chēng)為向下流動(dòng)的豎直槽道, 當(dāng)流動(dòng)方向與重力相反, 則稱(chēng)為向上流動(dòng)的豎直槽道, 如圖17. Challabotla 等(2016b, 2016c)和Yuan 等(2017)的研究表明重力會(huì)影響到顆粒的聚集形態(tài), 使得顆粒有顯著的法向輸運(yùn). 當(dāng)流動(dòng)方向向上時(shí), 非球形顆粒的分布會(huì)更加均勻, 而流動(dòng)方向向下時(shí), 顆粒聚集十分明顯, 而該影響隨著顆粒慣性的增大而減弱,且與顆粒形狀的關(guān)系不大. 除了非常靠近壁面的位置, 重力并不會(huì)對(duì)顆粒的取向與轉(zhuǎn)動(dòng)行為有明顯的影響. Yuan 等(2018)進(jìn)一步發(fā)現(xiàn)重力的作用會(huì)影響湍泳聚集(Turbophoresis)的機(jī)制且使得遠(yuǎn)離壁面的運(yùn)動(dòng)顆粒比靠近壁面的運(yùn)動(dòng)顆粒更聚集在流體的低速條帶. 如果忽略顆粒慣性的作用, Qiu 等(2019)發(fā)現(xiàn)的非球形顆粒在向上流動(dòng)的豎直槽道中傾向聚集在壁面, 而在向下流動(dòng)的豎直槽道中傾向聚集在槽道中部. 這說(shuō)明, 如果重力作用與基于湍泳機(jī)制聚集的作用相反, 那么顆粒的分布將會(huì)趨于均勻 (即向下流動(dòng)情況) ; 如果兩種作用相同, 那么會(huì)增強(qiáng)顆粒在壁面的聚集(即向上流動(dòng)情況). 而重力方向?qū)е虏煌姆ㄏ蜉斶\(yùn)特點(diǎn)可能與壁面附近的流動(dòng)上拋與下掃事件以及流動(dòng)結(jié)構(gòu)相關(guān). 其具體的影響值得進(jìn)一步深入研究.

6.2 有限尺寸顆粒

Do-Quang 等(2014) 對(duì)慣性較小的纖維使用全分辨方法研究其在槽道湍流中的行為, 并發(fā)現(xiàn)纖維會(huì)因?yàn)槠渑c壁面的碰撞效應(yīng)而傾向性的聚集在高速條帶, 且較長(zhǎng)的纖維在槽道中部有朝著展向的趨勢(shì), 隨著逐漸靠近壁面, 較長(zhǎng)的纖維在逐漸在剪切平面內(nèi)轉(zhuǎn)動(dòng). 而在壁面附近纖維主要朝著流向方向. Eshghinejadfard 等(2017, 2018)以及Zhu 等(2018)進(jìn)一步研究了中性懸浮桿狀與碟狀顆粒在槽道中的行為以及其與流體的相互作用. 他們發(fā)現(xiàn)在緩沖區(qū)時(shí), 顆粒的運(yùn)動(dòng)速度要快于流體速度, 且桿狀顆粒傾向性朝著流向而碟狀顆粒傾向性地朝著壁面法向. 而在遠(yuǎn)離壁面的位置, 顆粒的取向又趨于各向同性, 這與點(diǎn)顆粒方法得到的結(jié)論類(lèi)似. Eshghinejadfard 等(2019)通過(guò)改變顆粒與流體的密度比來(lái)改變顆粒的慣性, 其研究結(jié)果表明顆粒與流體間的密度比增大會(huì)影響流體與顆粒的統(tǒng)計(jì)行為. 其中顆粒的取向統(tǒng)計(jì)行為與點(diǎn)顆粒類(lèi)似, 隨著顆粒慣性的增大(密度比增加), 桿狀顆粒與流向、碟狀顆粒與法向的取向行為被削弱. Zhu 等(2020)研究軸對(duì)稱(chēng)桿狀顆粒與碟狀顆粒在豎直向上流動(dòng)的槽道湍流中的沉降問(wèn)題. 其結(jié)果顯示沉降的顆粒傾向于朝槽道中部遷移, 這種效應(yīng)對(duì)于球形顆粒最明顯. 而非球形顆粒在近壁附近依然使其最長(zhǎng)軸沿流向方向, 但在槽道中部由于沉降作用呈現(xiàn)出垂直流向的關(guān)系. 同時(shí), Zhu 等(2020)還發(fā)現(xiàn)顆粒傾向性會(huì)聚集在高速條帶.

7 顆粒取向與轉(zhuǎn)動(dòng)行為的力學(xué)機(jī)理

非球形顆粒盡管外形各向異性, 但是對(duì)于大部分問(wèn)題中顆粒的聚集行為與顆粒形狀關(guān)系不明顯, 聚集的機(jī)理基本可以由湍流泳動(dòng)的機(jī)制以及近壁相干結(jié)構(gòu)進(jìn)行解釋. 本文將主要討論由于形狀各向異性導(dǎo)致的復(fù)雜顆粒取向與轉(zhuǎn)動(dòng)行為背后的機(jī)理.

圖17豎直槽道示意圖.(a)向下流動(dòng), (b)向上流動(dòng)

7.1 渦量與變形率張量的影響

目前, 一些工作嘗試解釋各向同性湍流和壁湍流中的顆粒復(fù)雜行為的力學(xué)機(jī)理, 其中有關(guān)顆粒行為與渦量以及變形率張量的第二主軸的研究較為普遍(Pumir & Wilkinson 2011, Parsa et al. 2012, Byron et al. 2015). 比如, 在均勻各向同性湍流中, 無(wú)慣性的桿狀顆粒與流體的渦量方向具有較強(qiáng)的一致, 而碟狀顆粒垂直于渦量方向(Pumir & Wilkinson 2011, Byron et al. 2015).Pumir 和Wilkinson(2011)認(rèn)為桿狀顆粒之所以與渦量方向一樣, 是因?yàn)闇u量方向的演化方程與桿狀顆粒的運(yùn)動(dòng)方程具有類(lèi)似性, 其中渦量方向的方程只比桿狀顆粒取向方程多了額外的黏性擴(kuò)散項(xiàng). 因?yàn)闂U狀顆粒與渦量方向的關(guān)系, 在統(tǒng)計(jì)意義上桿狀顆粒的繞回轉(zhuǎn)軸轉(zhuǎn)動(dòng)的速率(Spinning)要強(qiáng)于碟狀顆粒, 但桿狀顆粒翻轉(zhuǎn)速率(Tumbling)要弱于碟狀顆粒(Parsa et al. 2012,Marcus et al. 2014, Gustavsson et al. 2014, Byron et al. 2015). 不過(guò), 直覺(jué)認(rèn)為桿狀顆粒應(yīng)該更加傾向性地朝著當(dāng)?shù)氐淖冃温蕪埩康淖畲罄旆较? 然而, 學(xué)者們(Pumir & Wilkinson 2011, Gustavsson et al. 2014, Guala et al. 2005, Chevillard et al. 2013, Ni et al. 2014)卻發(fā)現(xiàn)桿狀顆粒和流體渦量并不會(huì)朝著變形率張量的最大主軸方向. 相反, 桿狀顆粒與流體渦量與變形率張量的第二特征向量相關(guān)較強(qiáng). Xu 等(2011)發(fā)現(xiàn)并提出渦量會(huì)跟隨變形率張量的最大拉伸方向變化, 但是當(dāng)渦量轉(zhuǎn)動(dòng)到相應(yīng)位置時(shí), 原來(lái)的流體變形率張量的主軸方向已經(jīng)隨流體轉(zhuǎn)到其他方向. 在壁湍流中, 壁面附近的強(qiáng)剪切的使得壁面處的渦量方向一直朝著壁面展向方向, 但是無(wú)慣性桿狀顆粒與碟狀顆粒的回轉(zhuǎn)軸方向卻分別傾向性地朝著流向與壁面法向方向, 與渦量方向的夾角均接近直角, 如圖18(a)所示. 如果考慮顆粒慣性, 壁面的強(qiáng)渦量給顆粒提供了持續(xù)不斷的能量, 顆粒為了維持穩(wěn)定, 在自身轉(zhuǎn)動(dòng)慣量的作用下, 桿狀顆粒與碟狀顆粒的轉(zhuǎn)動(dòng)均傾向性地在流向與法向平面轉(zhuǎn)動(dòng)(Zhao et al. 2015, 2019). 此時(shí), 在顆粒慣性的作用下, 桿狀顆粒垂直于渦量但碟狀顆?；剞D(zhuǎn)軸卻與渦量對(duì)齊(Zhao et al. 2015, 2019), 如圖18(b). 但是, 因?yàn)樵诓鄣乐胁坎鄣乐胁康耐牧髭吔诰鶆蚋飨蛲缘耐牧鞯牧鲃?dòng)特征(Andersson et al. 2015), 顆粒的取向規(guī)律基本與各向同性湍流中的相同(Zhao et al. 2015, 2019).

圖18不同形狀和顆粒慣性的非球形顆?；剞D(zhuǎn)軸方向與渦量的夾角分布(Zhao et al. 2015). (a)槽道中部, (b)近壁區(qū)

由此得知, 顆粒取向行為與渦量和變形率的統(tǒng)計(jì)關(guān)系并不具備普適性, 比如, 桿狀顆粒傾向性朝著渦量方向僅在槽道中部或者均勻各向同性湍流中成立, 但在各向異性的近壁湍流附近桿狀顆粒卻與渦量方向垂直(Zhao et al. 2015). 而且, 顆粒在二維湍流或準(zhǔn)二維流動(dòng)中均不存在上述關(guān)系, 因?yàn)樵诙S流動(dòng)中渦量始終垂直于二維平面.

7.2 拉格朗日結(jié)構(gòu)

無(wú)論是分析顆粒取向與流體渦量以及變形率張量的關(guān)系, 還是顆粒取向與標(biāo)量梯度存在的聯(lián)系, 這些研究均考慮顆粒取向與歐拉觀(guān)點(diǎn)下的流場(chǎng)物理量進(jìn)行的比較分析. 實(shí)際上, 對(duì)于顆粒而言, 其方程均在拉格朗日觀(guān)點(diǎn)下進(jìn)行點(diǎn)顆粒追蹤求解. 對(duì)于無(wú)慣性的微小橢球顆粒, 顆粒軌跡就是流體跡線(xiàn), 而顆粒的轉(zhuǎn)動(dòng)可以由當(dāng)?shù)氐乃俣忍荻刃畔⑴cJeffery 方程直接得到. 如果僅采用顆粒與歐拉物理量(比如渦量與變形率張量) 進(jìn)行比較, 將得到的是顆粒與瞬時(shí)的歐拉流場(chǎng)的關(guān)系, 其并不能真實(shí)反映顆粒的拉格朗日行為特性.

如何從拉格朗日觀(guān)點(diǎn)研究顆粒與流體拉格朗日量進(jìn)行比較？在拉格朗日觀(guān)點(diǎn)下, 一個(gè)初始是球形的流體微團(tuán)沿著流體軌跡線(xiàn)會(huì)發(fā)生拉伸變形, 其中最強(qiáng)拉伸與壓縮方向可以被認(rèn)為是流體拉格朗日的拉伸與壓縮方向(Ni et al. 2014, Zhao & Andersson 2016). Parsa 等(2011)則通過(guò)實(shí)驗(yàn)數(shù)據(jù)獲得桿狀顆粒與流體拉格朗日拉伸方向之間存在明顯關(guān)系, 并展示了流體拉格朗日拉伸結(jié)構(gòu)對(duì)桿狀顆粒取向行為的影響 (如圖19(a)(b)).

圖19纖維與拉格朗日拉伸結(jié)構(gòu)的關(guān)系. (a)周期流動(dòng)(Parsa et al. 2011), (b)非周期流動(dòng)(Parsa et al.2011), (c)拉格朗日拉伸(黑色箭頭)與壓縮(紅色箭頭)與拉格朗日結(jié)構(gòu)的關(guān)系(Cui Z & Zhao 2021)

流體拉格朗日拉伸與桿狀顆粒取向存在怎樣的關(guān)系呢? 如果考慮一個(gè)細(xì)長(zhǎng)桿狀顆粒, 其可以被近似地認(rèn)為是一個(gè)物質(zhì)線(xiàn)段, 其方向與流體拉格朗日拉伸方向漸近一致(Pumir et al. 2011,Parsa et al. 2012, Batchelor 1952, Johnson et al. 2017). 在早期研究工作就發(fā)現(xiàn)物質(zhì)線(xiàn)段方向以及流體拉格朗日的拉伸方向就與渦量以及變形率張量的第二特征值方向有較強(qiáng)的關(guān)系(Girimaji& Pope 1990, Guala et al. 2005). 桿狀顆粒與拉格朗日拉伸方向以及物質(zhì)線(xiàn)段存在一些差別, 比如桿狀顆粒始終是有限長(zhǎng)寬比, 且不會(huì)在流體作用下發(fā)生拉伸與壓縮的變形. Ni 等(2014)驗(yàn)證了桿狀顆粒(長(zhǎng)寬比大于10)的取向可以用流體拉格朗日拉伸方向近似替代. 與此同時(shí), Zhao 和Andersson(2016)的研究結(jié)果也發(fā)現(xiàn)在槽道湍流中顆粒的朝向與Ni 等(2014)在各向同性湍流中的結(jié)論類(lèi)似, 即細(xì)長(zhǎng)桿狀顆粒的取向與流體拉格朗日最大拉伸方向有關(guān), 而扁平碟狀顆粒取向與流體拉格朗日最大壓縮方向相關(guān). 但是, Cui Z 等(2020)的研究表明顆粒形狀(長(zhǎng)寬比大于10)依然會(huì)對(duì)槽道的近壁附近的桿狀顆粒取向行為產(chǎn)生顯著影響. 而產(chǎn)生這一現(xiàn)象的原因可以歸結(jié)為槽道不同位置處流場(chǎng)的平均剪切與速度梯度脈動(dòng)強(qiáng)度的比值較大. 如果考慮桿狀顆粒在流體拉格朗日拉伸與壓縮方向組成的坐標(biāo)系中的取向分布, Cui Z 等(2020)的工作表明無(wú)論在壁面附近還是在槽道中部不同形狀的桿狀顆粒的取向分布均存在一個(gè)平臺(tái)區(qū)和冪律接近-2 的尾部區(qū).平臺(tái)區(qū)的存在就意味著顆粒在該區(qū)間內(nèi)取向分布相對(duì)均勻, 對(duì)應(yīng)著顆粒取向與流體拉格朗日拉伸方向的差別. 越靠近壁面, 平臺(tái)區(qū)的范圍更大, 顆粒取向與流體拉格朗日拉伸方向的差別就更明顯.類(lèi)似的平臺(tái)與冪律區(qū)的分布也意味著桿狀顆粒與流體拉格朗日拉伸方向的關(guān)系在全場(chǎng)是類(lèi)似的. 不過(guò), 越靠近壁面, 流場(chǎng)的平均剪切與速度梯度脈動(dòng)的比值變大, 分布函數(shù)中的平臺(tái)區(qū)變寬, 桿狀顆粒與流體拉格朗日拉伸方向的差別就變得更加顯著. Ni 等(2014)、Zhao 和Andersson (2016)以及Cui Z 等(2020)的研究均說(shuō)明了一個(gè)共同的規(guī)律, 即隨著長(zhǎng)寬比的增加/趨于零, 桿狀/碟狀顆粒逐漸與流體拉格朗日拉伸/壓縮方向一致. 只不過(guò)在剪切湍流中, 若想讓顆粒與流體拉格朗日方向完全一致, 則需要更小的形狀差異δΛ=1-|Λ| . 有關(guān)這一部分的闡述可見(jiàn)崔智文和趙立豪(2021)關(guān)于該內(nèi)容的詳細(xì)綜述.

隨著人們對(duì)無(wú)慣性橢球顆粒取向行為研究的深入, 流體的拉格朗日拉伸與壓縮方向與非球形顆粒取向之間的數(shù)學(xué)關(guān)系也逐漸明晰. 在Pumir 等(2011)、Ni 等(2014)以及Cui Z 等(2020,2021)的工作中均將Jeffery 方程中的Λ=1 或 - 1 去近似替代流體拉格朗日拉伸或壓縮以及細(xì)長(zhǎng)的桿狀或扁平的碟狀顆粒. 其背后的數(shù)學(xué)關(guān)系可以由Batchelor (1952)無(wú)限小物質(zhì)線(xiàn)與物質(zhì)面的演化得到. 不過(guò) Balkovsky 和Fouxon (1999)的工作發(fā)現(xiàn)當(dāng)積分時(shí)間足夠長(zhǎng)后左柯西格林張量的最大與最小的特征向量(即流體拉格朗日拉伸與壓縮方向)的方程形式與Jeffery 方程中令Λ=1或 - 1 一樣. Cui Z 和Zhao (2021)則在該工作的基礎(chǔ)上, 進(jìn)一步推導(dǎo)并完善了左柯西格林張量特征方向的時(shí)間演化方程, 并發(fā)現(xiàn)其在數(shù)學(xué)形式上與無(wú)慣性三軸不等橢球顆粒主軸的時(shí)間演化方程的高度相似性, 見(jiàn)式(49)和式(50). 在公式中, 隨著積分時(shí)間的延長(zhǎng), 柯西格林張量方程中的參數(shù)?i(i=1,2,3) 將 由無(wú)窮趨于 1,-1 和 1 . 這恰好與長(zhǎng)軸比中間軸無(wú)限大、中間軸比短軸無(wú)限大的三軸不等橢球顆粒的取向方程一致, 即Λ1=1,Λ2=-1,Λ3=1 . 此時(shí), 式(49)與式(50)完全一致. 隨后Cui Z 和Zhao (2021)利用該特性結(jié)合三軸不等橢球顆粒求解過(guò)程提出了一種新穎的左柯西格林張量的方法, 詳細(xì)內(nèi)容可參見(jiàn)(Cui Z & Zhao 2021).

左柯西格林張量三個(gè)特征向量〈e1,e2,e3〉的時(shí)間演化方程

對(duì)應(yīng)的三軸不等橢球顆粒的主軸〈p1,p2,p3〉的時(shí)間演化方程

8 顆粒在湍流中的減阻控制作用

湍流的減阻控制在工業(yè)生產(chǎn)和交通輸運(yùn)等工程領(lǐng)域具有重要意義. 大量實(shí)驗(yàn)發(fā)現(xiàn), 湍流里放入少量的添加劑, 如柔性聚合物、表面活化劑、氣泡等, 可能產(chǎn)生一定的減阻效應(yīng). 其中柔性聚合物減阻明顯, 且使用較為方便, 在許多領(lǐng)域得到了廣泛應(yīng)用. 但其不受強(qiáng)剪切, 不適用于剪切較強(qiáng)的流動(dòng)中(Reddy & Singh 1985). 部分研究者通過(guò)實(shí)驗(yàn)發(fā)現(xiàn)剛性纖維可經(jīng)受強(qiáng)剪切, 且同樣具有減阻作用(Radin et al. 1975, Gyr & Bewersdorff 1995). 在數(shù)值模擬和理論分析中, 細(xì)小的剛性纖維通常被近似為無(wú)慣性的細(xì)長(zhǎng)橢球顆粒 (如Brenner 1974, Gillissen et al. 2008) . 無(wú)慣性橢球顆粒通過(guò)額外的應(yīng)力影響附近流體, 而反饋力和力矩則可忽略(Guazzelli & Morris 2011). Den Toonder 等(1997)首先利用直接數(shù)值模擬方法研究含細(xì)長(zhǎng)顆粒的圓管湍流, 盡管他們簡(jiǎn)單地假定顆粒的取向沿著當(dāng)?shù)氐牧黧w速度向量的方向, 但仍然發(fā)現(xiàn)細(xì)長(zhǎng)的剛性顆粒具有湍流減阻效應(yīng),且湍流統(tǒng)計(jì)結(jié)果與柔性聚合物的實(shí)驗(yàn)結(jié)果在定性上基本一致. Paschkewitz 等(2004)通過(guò)歐拉-歐拉雙向耦合方法研究了最小尺寸槽道湍流中細(xì)長(zhǎng)顆粒的減阻作用, 探討了細(xì)長(zhǎng)顆粒的布朗運(yùn)動(dòng)、體積濃度、長(zhǎng)寬比和顆粒取向的封閉模型對(duì)減阻強(qiáng)弱的影響, 實(shí)現(xiàn)最大減阻率可達(dá)26%, 還提出了一種可能的湍流減阻機(jī)理. 他們認(rèn)為細(xì)長(zhǎng)顆粒在準(zhǔn)流向渦之間的拉伸區(qū)會(huì)產(chǎn)生強(qiáng)烈的應(yīng)力脈動(dòng), 從而減弱了近壁渦旋結(jié)構(gòu)和改變了近壁湍流的自維持過(guò)程, 最終導(dǎo)致了湍流減阻現(xiàn)象,如圖20 所示. Gillissen 等(2008)進(jìn)一步發(fā)現(xiàn)細(xì)長(zhǎng)顆粒在不同流動(dòng)雷諾數(shù)下都具有湍流減阻效應(yīng).在這些數(shù)值模擬中, 顆粒相是在歐拉觀(guān)點(diǎn)下結(jié)合封閉模型進(jìn)行求解, 而非使用拉格朗日點(diǎn)顆粒法直接求解單個(gè)顆粒的動(dòng)力學(xué)方程. Moosaie 和Manhart (2013)首次通過(guò)歐拉-拉格朗日雙向耦合方法數(shù)值模擬了無(wú)慣性細(xì)長(zhǎng)橢球顆粒與槽道湍流的相互作用, 同樣發(fā)現(xiàn)了壁湍流減阻的現(xiàn)象和流動(dòng)特征. 不過(guò), 相比于細(xì)長(zhǎng)顆粒, 扁平顆粒對(duì)壁湍流的影響的相關(guān)研究則非常少. 單向耦合的數(shù)值模擬結(jié)果表明, 長(zhǎng)寬比較小的扁平顆粒也會(huì)產(chǎn)生強(qiáng)烈的顆粒應(yīng)力脈動(dòng), 這意味著一定濃度的扁平顆?；蛟S能夠調(diào)制湍流(Wang & Zhao 2020). Wang 等(2021)最近利用雙向耦合的直接數(shù)值模擬, 研究了體積分?jǐn)?shù)為0.75%的無(wú)慣性扁平顆粒和細(xì)長(zhǎng)顆粒對(duì)槽道湍流的影響, 發(fā)現(xiàn)長(zhǎng)寬比為100 的細(xì)長(zhǎng)顆粒減阻14.93%, 長(zhǎng)寬比為0.01 和0.002 的扁平顆粒分別減阻1.92%和7.15%. 他們提出了無(wú)慣性軸對(duì)稱(chēng)橢球顆粒的湍流減阻機(jī)理, 認(rèn)為顆粒通過(guò)在展向和壁面法向?qū)α黧w做的功抑制了準(zhǔn)流向渦結(jié)構(gòu), 從而減弱了上拋下掃事件的強(qiáng)度, 導(dǎo)致了雷諾切應(yīng)力的減小, 同時(shí)顆粒切應(yīng)力會(huì)部分補(bǔ)償雷諾切應(yīng)力的減小, 兩者之和的加權(quán)積分確定了阻力系數(shù)的大小. 對(duì)于扁平顆粒而言, 雖然雷諾切應(yīng)力出現(xiàn)了明顯減小, 但近壁區(qū)處較大的顆粒切應(yīng)力減弱了湍流減阻的效果.

另一方面, 一些學(xué)者利用全分辨方法數(shù)值模擬含有有限尺寸顆粒的壁湍流, 發(fā)現(xiàn)有限尺寸的橢球顆粒也可能減小流動(dòng)阻力. Ardekani 等(2017b)發(fā)現(xiàn)中性懸浮的扁平顆粒在近壁區(qū)轉(zhuǎn)動(dòng)較慢且回轉(zhuǎn)軸傾向沿壁面法向, 這導(dǎo)致了流體雷諾應(yīng)力和湍流產(chǎn)生項(xiàng)的減小, 從而產(chǎn)生湍流減阻.Eshghinejadfard 等(2018)觀(guān)察到體積分?jǐn)?shù)為10%的有限尺寸的扁平顆粒使得流體平均速度增加了1.3%. 他們提出, 展向速度脈動(dòng)的減小和條帶結(jié)構(gòu)的展向間距增大是扁平顆粒減阻的主要原因, 而顆粒對(duì)渦結(jié)構(gòu)的增強(qiáng)作用則抑制了減阻效果. Zhu 等(2018)數(shù)值模擬了含有限尺寸的細(xì)長(zhǎng)和扁平顆粒的槽道湍流, 發(fā)現(xiàn)顆粒偏離球形的程度越大減阻越明顯, 他們認(rèn)為近壁區(qū)較低的顆粒體積分?jǐn)?shù)和較小的流體雷諾應(yīng)力導(dǎo)致了湍流減阻.

9 總結(jié)與展望

本文主要綜述了非球形顆粒兩相流的數(shù)值模擬方法、非球形顆粒在不同類(lèi)型的簡(jiǎn)單流場(chǎng)與復(fù)雜湍流場(chǎng)中的運(yùn)動(dòng)行為的研究進(jìn)展, 并對(duì)非球形顆粒取向與轉(zhuǎn)動(dòng)行為背后的機(jī)理進(jìn)行了闡述與討論. 根據(jù)目前的研究現(xiàn)狀, 未來(lái)非球形顆粒兩相流可能的研究方向包括:

圖20纖維減阻機(jī)制示意圖(Paschkewitz et al. 2004)

(1) 非球形顆粒在流體中受力模型的建立與完善. 在目前主流的非球形顆粒兩相流的點(diǎn)顆粒模型中, 顆粒的阻力模型、力矩模型以及顆粒對(duì)流體反饋?zhàn)饔媚Ｐ椭饕獊?lái)自20 世紀(jì)Jeffery,Batchelor, Brenner 以及他們合作者的相關(guān)系列工作. 近年Dabade 等(2015)的流體慣性力矩模型是對(duì)非球形顆粒理論模型的重要補(bǔ)充, 但該模型的正確性與普適性還有待數(shù)值模擬以及實(shí)驗(yàn)的驗(yàn)證. 不過(guò)與球形顆粒相比, 非球形顆粒的模型依然還不夠完善, 比如Maxey-Riley 模型中提及的壓力梯度力、附加質(zhì)量力以及Basset 歷史力等模型在非球形顆粒問(wèn)題中往往被人為地忽略. 對(duì)于非球形顆粒而言, 這些力的準(zhǔn)確表達(dá)式以及其對(duì)非球形顆粒運(yùn)動(dòng)的具體影響尚不明確.目前, 除了從理論層面推導(dǎo), 未來(lái)可主要通過(guò)全分辨數(shù)值模擬的手段研究單顆粒在簡(jiǎn)單流場(chǎng)中的運(yùn)動(dòng)來(lái)驗(yàn)證與分析相關(guān)模型的正確性, 或進(jìn)行相關(guān)模型的建立與改進(jìn)工作. 非球形顆粒模型的建立與完善, 一方面能幫助研究者對(duì)非球形顆粒行為進(jìn)行精準(zhǔn)捕捉, 進(jìn)而改進(jìn)目前的顆粒模擬方法, 從而進(jìn)行大規(guī)模的非球形顆粒兩相流的數(shù)值模擬, 另一方面能幫助研究者解釋非球形顆粒復(fù)雜行為背后的機(jī)理, 為相關(guān)工業(yè)設(shè)計(jì)提供理論依據(jù).

(2) 中高雷諾數(shù)非球形顆粒湍流兩相流的研究. 目前絕大部分非球形顆粒兩相流的數(shù)值模擬主要集中在較低雷諾數(shù)的情況, 有關(guān)非球形顆粒在中高等雷諾數(shù)的情況的研究屈指可數(shù), 以非球形顆粒槽道湍流為例, 目前, 摩擦雷諾數(shù)Reτ超過(guò)1000 的直接數(shù)值模擬的工作僅Jie 等(2019)與Ouchene 等(2018)兩項(xiàng). 然而在實(shí)際的科學(xué)和工業(yè)問(wèn)題中, 流場(chǎng)往往具有較高的雷諾數(shù), 在高雷諾數(shù)下, 流動(dòng)中會(huì)存在明顯的大尺度結(jié)構(gòu). 大尺度結(jié)構(gòu)對(duì)非球形顆粒的作用依然有待進(jìn)一步研究. 與此同時(shí), 對(duì)于高雷諾數(shù)的流場(chǎng)研究已經(jīng)有比較成熟的大渦模擬的方法取代直接數(shù)值模擬.但大渦模擬是對(duì)小尺度的流動(dòng)進(jìn)行模化, 而非球形顆粒的運(yùn)動(dòng)行為又往往對(duì)小尺度的流動(dòng)結(jié)構(gòu)比較敏感, Chen 等(2016)的工作展示了不同流體的亞格子模式就會(huì)對(duì)非球形顆粒的轉(zhuǎn)動(dòng)與取向行為產(chǎn)生顯著影響. 因此, 未來(lái)對(duì)非球形顆粒相關(guān)的亞格子模型以及對(duì)現(xiàn)有的湍流亞格子模型的改進(jìn)依然是一個(gè)比較重要且十分具有挑戰(zhàn)性的難題. 而最近的Pujara 等(2021)的工作表明無(wú)慣性橢球顆粒的運(yùn)動(dòng)行為具有形狀與尺度的無(wú)關(guān)性, 這可能展示了湍流的某種拉格朗日不變量隨湍流尺度無(wú)關(guān). 這或許可以為非球形顆粒大渦模擬的研究工作提供新的研究思路.

(3) 非球形顆粒運(yùn)動(dòng)行為的機(jī)理分析. 非球形顆粒在流動(dòng)中的復(fù)雜運(yùn)動(dòng)行為一方面由顆粒自身復(fù)雜的外形、顆粒慣性、顆粒運(yùn)動(dòng)學(xué)模型引起, 另一方面也受到顆粒周?chē)膹?fù)雜的流動(dòng)環(huán)境的影響. 目前有關(guān)非球形顆粒的機(jī)理研究存在兩種基本方式, 第一種就是用當(dāng)?shù)氐牧鲌?chǎng)信息與顆粒運(yùn)動(dòng)行為建立對(duì)應(yīng)關(guān)系, 比如渦量與變形率張量與顆粒行為的聯(lián)系. 該方法的最大優(yōu)點(diǎn)是直觀(guān), 但不可避免地忽略了顆粒拉格朗日運(yùn)動(dòng)的影響. 因此, 該類(lèi)方法的最大缺點(diǎn)在于缺少拉格朗日信息, 且利用渦量與變形率解釋顆粒取向并不具備普適性. 盡管該類(lèi)方法在各向同性湍流中取得非常大的成功, 但是我們還是能看到其在具有強(qiáng)剪切的壁面以及準(zhǔn)二維的流動(dòng)中是基本上失效的. 而第二種方式就是利用拉格朗日拉伸與壓縮方向來(lái)解釋顆粒的取向行為, 這是從拉格朗日的觀(guān)點(diǎn)出發(fā)進(jìn)行解釋, 且該方法包括了必要的拉格朗日信息, 因此橢球顆粒的取向與流體拉格朗日拉伸或壓縮方向表現(xiàn)出極強(qiáng)的相關(guān)性. 然而, 該方法通常計(jì)算復(fù)雜, 需要得到大量的流體拉格朗日軌跡線(xiàn), 且一般僅對(duì)于完全跟隨流體質(zhì)點(diǎn)的非球形顆粒有比較好的效果. 非球形顆粒運(yùn)動(dòng)行為的機(jī)理研究的主要目的之一就是能夠希望通過(guò)周?chē)牧鲌?chǎng)信息推測(cè)出非球形顆粒的可能運(yùn)動(dòng)行為, 從而可以實(shí)現(xiàn)顆粒行為的控制或者預(yù)測(cè). 然而, 在實(shí)際問(wèn)題中流場(chǎng)的情況非常復(fù)雜, 使用第二種方式分析和預(yù)測(cè)顆粒的取向行為需要花費(fèi)較大的計(jì)算成本. 因此, 第一種方式將是比較經(jīng)濟(jì)且有效的. 那么如何在第一種方式的基礎(chǔ)上或者其他方式, 提出更加有效的分析方法將可能是未來(lái)非球形顆粒兩相流研究一個(gè)重要的方向.

(4) 復(fù)雜流場(chǎng)中復(fù)雜顆粒的運(yùn)動(dòng)行為. 本綜述主要介紹了非球形顆粒在不可壓牛頓流體中的運(yùn)動(dòng)行為, 且流動(dòng)的幾何邊界非常簡(jiǎn)單. 未來(lái)可以從兩個(gè)方面來(lái)拓展非球形顆粒兩相流的研究范圍. 一個(gè)方面, 可以通過(guò)改變流場(chǎng)屬性、幾何等特征, 比如非牛頓流、密度分層流、熱對(duì)流或者可壓縮流動(dòng), 或研究粗糙壁面湍流、射流等復(fù)雜流動(dòng)問(wèn)題中的顆粒行為. 另一方面, 可以通過(guò)改變顆粒模型進(jìn)而研究不同非球形顆粒模型對(duì)顆粒運(yùn)動(dòng)行為的影響, 比如游動(dòng)、顆粒形狀反對(duì)稱(chēng)、質(zhì)心偏置等. 另外, 有關(guān)結(jié)合機(jī)器學(xué)習(xí)的智能顆粒的研究在近些年也得到了廣泛地關(guān)注, 相關(guān)綜述可見(jiàn)(邱敬然和趙立豪 2021).

致謝

國(guó)家自然科學(xué)基金(11 702 158, 11 911 530 141)和清華大學(xué)國(guó)強(qiáng)研究院(2019GQG1012)資助項(xiàng)目.