微納顆粒兩相流系統中湍流脈動效應研究綜述

楊海林，林建忠

(浙江大學 航空航天學院 流體工程研究所，杭州 310027)

0 引 言

在自然界和工程實際應用的氣固兩相流中，固相的顆粒如粉、塵、煙、霾等不僅受到流體相的輸運作用，其自身還往往出現由物理或化學作用導致的生長、凝并（顆粒碰撞后形成一個顆粒）、破碎等現象。顆粒一般動力學方程（General dynamic equation, GDE），又稱顆粒數平衡方程（Population balance equation, PBE），可用來描述顆粒相數密度在流體輸運下出現以上現象后的演變過程[1]。

在大多數情況下，兩相流處于湍流狀態，湍流導致的流動無序性與多尺度特性，使得兩相湍流變得很復雜，因而更具有研究價值。流場物理量的脈動是湍流場的一個重要特征，對于湍流脈動特性的研究，其科學意義不僅局限于湍流或多相流領域[2]。本文關注微納顆粒兩相湍流中的湍流脈動效應，在微納顆粒兩相湍流中，顆粒的形成與生長受到環境因素如流場速度、流場剪切率、化學組分濃度、溫度、壓力等的影響，而這些環境因素會受到湍流脈動的影響變化，所以顆粒的形成和生長速率會受到湍流脈動的影響。于是，若不考慮湍流的脈動而僅僅用平均后的環境因素來計算顆粒的形成和生長就會與實際情形產生偏差。例如，在云層凝結的過程中，湍流導致的溫度脈動會使水蒸氣飽和度發生變化，從而使得用平均參數給出的云霧顆粒形成率與生長率與實際不符[3-8]。在實驗室條件下，湍流造成的溫度脈動、物質組分脈動也不可忽略[9-11]，否則會使實驗測量結果產生偏差[12-13]。在數值模擬方面，已有一些探討湍流脈動影響的研究，例如將DNS得到的結果與用平均場求得的結果進行對比[14-18]，將LES得到的結果與用雷諾平均方法得到的結果進行對比[19-24]，將用和不用脈動量概率密度函數進行修正的結果進行對比[25-27]。除了顆粒的形成與生長受湍流脈動的影響外，顆粒的凝并與破碎也受湍流脈動的影響，因為顆粒的布朗凝并、湍流凝并、湍流破碎都與顆粒的濃度有關，而顆粒濃度在湍流場中也存在脈動[2,10,28-34]。湍流脈動對顆粒凝并影響的研究雖然較少，但這種影響是存在的[25,35-38]。

1 顆粒相方程

1.1 顆粒一般動力學方程

顆粒一般動力學方程是描述顆粒數密度函數（Number density function, NDF）的控制方程。假設顆粒尺度小于最小流動尺度且顆粒雷諾數遠小于1，顆粒數密度分布n(ξ;x,t)表 示顆粒在不同尺度ξ 、時間t與空間x上的數量分布。根據需要，ξ可以選取顆粒體積v或顆粒直徑dp，本文選顆粒體積，簡寫作n(v)或n，則一般動力學方程為[1]：

式中u為 流場速度，Dp為顆粒擴散系數，方程左邊依次為非定常項、對流項和擴散項，方程右邊表示導致顆粒數密度變化的因素，依次為成核項、凝積生長項（小顆粒沉積在大顆粒上使大顆粒尺度增大）、碰撞凝并項和破碎項。方程（1）是通過對顆粒數密度建立守恒方程而得到，具有合理性，適用于數密度連續且可微的情形。

（1）成核項。成核是指通過結晶或其他化學作用形成初始體積為v0的顆粒的過程：

式中B(Y) 為顆粒生成核函數，Y為各類環境參數的向量，根據顆粒的性質不同，環境參數向量包含的變量也不一樣，例如水蒸氣中液滴核的生成核函數中包括溫度、壓力、表面張力和飽和度等變量，而一些化學變化中則包含各類化學組分的變量。

（2）生長項。顆粒的凝積生長是指其他物質凝結、沉積在顆粒上，使顆粒尺度增長的過程：

式中G(Y,v)為顆粒生長核函數，也具有多種形式。此外，一些顆粒尺度減小的過程如蒸發、升華，也可用生長項描述，只是其值為負數。

（3）凝并項。凝并是指顆粒間相互碰撞聚集導致顆粒整體尺度增大、數密度減少的過程：

式中 β(v,v1)為 凝并核函數，描述體積為v和v1兩種顆粒的碰撞概率。

（4）破碎項。破碎是指顆粒因無法支撐團聚狀態而導致裂解，使得顆粒整體尺度減小、數密度增加的過程：

式中a(v)是 破碎核函數，描述體積為v的顆粒的破碎概率，b(v|v1) 是破碎子函數，描述體積為v1的顆粒在發生破碎時其子顆粒體積為v的條件概率。

對方程（1）的求解難點在于對源項的處理，當對每個源項給出相應的表示方式后，便可對方程（1）進行求解。除了極少數情況能得到方程（1）的相似解[1]，大多數情況需要采用數值方法求解。在數值求解方法中，分區法和節點法是在將顆粒尺度離散化的基礎上得到離散的尺度分布，然后通過計算大量離散尺度下的顆粒數密度得到收斂的數密度分布[39-40]。矩方法則是將顆粒數密度分布對顆粒尺度 ξk進行積分，經過近似假設和數學處理得到k階矩方程后進行求解，該方法大大簡化了求解的方程，且不同階的矩具有特定的物理意義，因而被大量使用[41-46]。此外，還有運用統計學方法，將顆粒的尺度分布處理成聯合PDF方程進行求解[25-27,47-50]，也有對顆粒的虛擬拉格朗日粒子進行蒙特卡洛模擬的方法[51-55]。前者的優點是能利用比較成熟的微積分方程數值求解方法進行求解，缺點是聯合概率密度函數的確定具有一定的復雜性。后者的優點是不必進行離散化處理，缺點是往往需要較大的計算量。

1.2 湍流場的顆粒一般動力學方程

采用雷諾平均方法，將流場速度和顆粒數密度表示成平均量和脈動量之和：

將其代入方程(1)后進行平均得：

左邊第三項為脈動對流項，通常采用梯度擴散假說進行封閉，引入湍流擴散系數Dt，使

將其與左邊第四項合并得：

則方程左邊完全封閉，接下來處理右邊。

(1) 成核項。經平均后的成核項為：

(Y)為平均成核率。有些研究采取簡化方案，通過平均環境參數Yˉ 求出表觀成核率B() 來代替(Y)，但通常

(2)生長項。生長項的雷諾平均有幾種情況：

①最簡單的情況是生長核函數與環境參數無關，用G(v)表示，則生長項的雷諾平均為：

該情況的生長項無需封閉，當溫度、化學組分等環境參數對顆粒生長的影響很小時，可直接使用。

②若生長核函數與環境參數相關，但與顆粒尺度無關，核函數會受到環境參數脈動的影響，則生長核函數可寫成G(Y)的 形式，對求平均后，通常生長核函數為環境變量的非線性函數，即因此有

式中的實際平均生長率Gˉ(Y)與表觀生長率不同，即，需要模化。當顆粒尺度的演變局限在某一范圍內且生長核函數可近似視為與顆粒尺度無關時，可用上式。

③若生長核函數與環境參數以及顆粒尺度v相關，則：

(3) 凝并項。該項的雷諾平均也分有幾種情況：

①凝并核函數 β(v,v1)與環境參數無關或對環境參數變化不敏感，則：

②凝并核函數 β(v,v1)與環境參數有關，則環境參數的脈動對 β(v,v1)有影響需 要模化，此時要根據具體情況對

進行展開，例如湍流凝并核函數為[56]:

式中α為常數，剪切率 Γη與 湍流耗散率 ε和流體黏度v相關，在Kolmogorov尺度上 Γη=(ε/v)1/2。

(4) 破碎項。

①當破碎核函數不含任何湍流脈動量時，平均后的破碎項為：

②當破碎核函數含湍流脈變量時，平均后的破碎項為：

如湍流破碎核為：

式中：kb為 玻爾茲曼常數；μ 為流體動力黏度； Γ*是特征切應力；指數q為正實數，由實驗確定。

湍流脈動引起顆粒的成核、生長、凝并與破碎項使方程(1)變得不封閉，目前尚無對這些項模化的公認模型，大部分數值模擬中只能直接忽略這些項，而已有研究表明，有些情況下這些項并不可忽略，例如湍流微尺度與顆粒尺度相當的情形[1]。

2 湍流脈動導致的顆粒成核與生長

隨著對湍流研究的不斷深入，湍流場中脈動壓力、脈動溫度對顆粒成核與生長的影響逐漸引起人們的關注。最早開始研究湍流脈動對顆粒生成與生長影響的是Levin和Sedunov[57]，他們認為湍流場中壓力的脈動影響飽和度的穩定，從而改變雨滴的生長率，使雨滴的尺度分布變寬。然而，他們沒有給出顆粒濃度和粒徑分布的最終信息。隨著對概率密度函數(PDF)方程求解技術的提高，有關湍流脈動對顆粒生成與生長影響的研究有了比較有效的手段。Lesniewski和Friedlander[12]計算了圓湍射流場蒸汽冷凝過程中顆粒成核的速率，并對比了用溫度PDF和平均溫度場計算得出的結果。經典理論表明，在過飽和點的位置附近，顆粒成核速率會急劇增加，而在該位置之外成核速率幾乎為零。采用溫度PDF研究的結果表明，所有情況下都存在蒸汽冷凝成核，只是最大成核速率沒有經典理論給出的那么大。在考慮溫度脈動的情況下，得到的顆粒濃度分布較均勻且更寬。后來，Lesniewski和他的同事[13,58]繼續用二丁基鄰苯二甲酸與氮氣在圓湍射流中進行顆粒成核實驗，以說明脈動溫度與脈動組分濃度對成核的影響，結果表明湍流導致的脈動對顆粒生成與生長的影響不可忽略，這一結果引起了人們的重視，后續的數值模擬研究[15-17,19-22,26-27,37]中都采用了與該實驗相同的條件。

2.1 顆粒成核率、生長率的數值模擬結果

采用不同的數值模擬方法，得到的顆粒成核率、生長率的結果也不同。以下采取統一的表達方式來描述瞬時、表觀和各類平均的數值結果。以經典的蒸汽凝結均質成核速率表達式為例：

式中Pv、Nv表示蒸汽壓力和濃度，vm為 摩爾體積，m為分子質量，k為玻爾茲曼常數，T為絕對溫度，S為飽和度，σ表示液滴的表面張力。將各種環境參數表示為向量Y，則Y在湍流場中存在脈動。如果將用DNS得到的瞬時環境參數代入，得到的是瞬時成核速率B(Y)。將DNS給出的實時成核速率進行時間平均或系綜平均，可得到真實平均值若將各類環境參數的平均值代入，得到的是表觀值

要考慮環境參數脈動對成核的影響，一種方法是引入脈動環境參數已知的PDF，然后計算PDF的平均值。例如，考慮溫度和蒸汽濃度的脈動時，PDF平均成核速率為[13,19]：

式中 pdf〈〉是 溫度T和蒸汽濃度Nv的概率平均值，由其概率分布函數PDF得出。

Fager等[20]將用LES得到的環境參數代入成核速率式，得到空間濾波的成核速率，還以前面的成核率為例：

對生長核函數G(Y,v)而言，可以用上述的方法定義表觀值平均值空間濾波值G(〈Y〉,v)和 真實平均值

Veroli和Rigopoulos[27]采用蒙特卡洛方法計算了二丁基鄰苯二甲酸蒸汽在圓湍射流中的凝結成核，并對用平均溫度場和溫度場PDF得到的表觀飽和度和 PDF平均飽和度S(pdf〈T〉)進行了比較，發現在射 流 出 口 附 近，S()的 值 更 高，分 布 更 加 集 中，而S(pdf〈T〉)的 值 較 低，分 布 較 為 分 散，這 一 結 果 與Lesniewski和Friedlander的實驗結果[13]相符。

Zhou等[17]用DNS方法模擬流場、積分矩方法(Quadrature Method of Moment, QMOM)求解顆粒相的源項，計算了二丁基鄰苯二甲酸蒸汽在混合剪切層中的冷凝成核與凝積生長，并將成核速率和生長速率的真 實 平 均 值與 表 觀 平 均 值進行了比較，結果如圖1所示。可見在考慮湍流脈動的情況下，成核速率在空間的分布更分散，最大值更小，而生長速率在空間的分布趨勢沒有明顯變化，只是總體數值降低。

圖1 成核速率與粒徑凝積生長速率在混合剪切層的空間分布[17]Fig. 1 Space distribution of nucleation rate and condensation growth rate in mixing shear layer[17]

Pesmazoglou等[22]用LES數值模擬了圓湍射流中二丁基鄰苯二甲酸蒸汽的凝結成核，并將由LES和RANS方法得到的空間濾波成核速率B(〈Y〉)和表觀成核速率進行了比較，得到了類似圖1（a）的結果，說明用兩種方法得到的結果存在差異。

2.2 湍流脈動對顆粒成核率、生長率的影響程度

在對受湍流脈動影響的顆粒成核項和生長項進行模化時，首先需確定湍流脈動在特定條件下對這兩項的影響程度。例如，在蒸汽凝結成核過程中，顆粒成核率對溫度的變化很敏感，Garrick等對二丁基鄰苯二甲酸蒸汽的凝結成核進行了數值模擬[15-16,20-21]，發現湍流場中的顆粒成核率主要受溫度、飽和度、組分濃度這三者脈動的影響。

圖2 圖2 濾波結果與瞬時值比較的散點圖[15]Fig. 2 Scatter plots of filtered results vs instantaneous results[15]

進一步對數據進行處理和研究，可以得到不同溫度脈動情況下的成核率亞格子尺度脈動，結果如圖3所示。對比圖3（a）和圖3（b）可知，成核率的負脈動比正脈動的幅值大，溫度的負脈動和正脈動對成核率的負脈動和正脈動的影響最大，成核率約為1 ×1012m3/s時對溫度脈動最敏感。

圖3 不同標準化脈動溫度θ ′條件下，SGS尺度脈動成核率的概率分布函數[15]Fig. 3 Under different normalized fluctuation temperaturecondition, the PDF of SGS nucleation rate [15]

3 湍流脈動導致的顆粒凝并與破碎

在顆粒凝并項方程(13)中，包含了n′(v)n′(v1)，該項是顆粒脈動數密度的二階關聯，是個未知項，需要模化[59-60]。由方程(13)可知，由湍流脈動引起的顆粒凝并與顆粒的數密度脈動有關。Becker等[28,61-62]由實驗證明了顆粒濃度在湍流場中確實存在脈動，并測量給出了顆粒濃度脈動均方根、脈動強度、脈動譜在圓管湍射流場中的分布。文獻[63]對相關的實驗研究進行了總結，在此不贅述。與雷諾應力的表達式有相似之處，只是后者為脈動速度，于是人們尋找與湍動能和湍流耗散的關系。Spalding[29]曾提出關于顆粒濃度脈動的預測方程，該方程將顆粒濃度脈動的均根方值與湍動能的生成與耗散相聯系，最終經過修正得到k-ε-g方程[64]。

隨著空間濾波方法的出現以及LES方法的發展，人們提出了顆粒數密度脈動的空間自關聯形式〈n′(x)n′(x+r)〉， 引入顆粒徑向分布函數g(r)將脈動項與平均項聯系起來，并對g(r)的性質進行了研究[30-31,65-67]。Eaton和Fessler[32]認為在小的時間尺度內，湍流渦的離心作用將顆粒甩到流場的高應變區域，形成了顆粒數密度聚集區，該區的范圍與顆粒Stokes數有關。Stokes數定義為顆粒松弛時間和流體特征時間之比：

式中 τp和 τη分別為顆粒松弛時間和流體特征時間，τη=(v/ε)1/2為Kolmogorov時 間 尺 度， η=(v3/ε)1/4為Kolmogorov空間尺度， ρp和 ρf分別為顆粒密度和流體密度。

3.1 顆粒聚集效應與凝并核函數的修正

為了體現湍流脈動對顆粒凝并的影響，一種思路是對顆粒凝并率進行湍流脈動的修正，即給出湍流凝并核函數β。基于這種思路，可以將顆粒徑向分布函數g(r)引入作為增強系數[67-73]，從而得到聚集效應下的凝并核函數。Reade和Collins[69]在DNS方法中引入體現顆粒數密度脈動對顆粒凝并影響的增強因子，并用函數G11和G12來分別描述單分散和雙分散顆粒的情形。Zhou等[71]則從歐拉方法和統計方法出發提出增強因子：

式 中N1和N2分 別 代 表 半 徑 為a1和a2顆 粒 的 數 密 度。通過DNS可計算出增強因子，從而得到適用于高顆粒數密度條件下的修正凝并核函數。

Sundaram[67]提出顆粒聚集時單分散顆粒的凝并核函數：

式中wr為 兩個顆粒的徑向相對速度，P(wr|r)為兩個顆粒在距離為r時其徑向相對速度為wr的條件概率。

基于上述凝并核函數的形式，Reade和Collins[69]得到顆粒徑向分布函數為St)， 該函數與粒徑、Stokes數相關，當St→0 或St→∞時，增強因子趨近于1，意味著湍流脈動導致的顆粒聚集效應不存在，只有在Stokes數達到1的量級時，增強因子達到最大。

Wang等[70]進一步得到新的凝并核函數β=這里的wr假設為高斯分布，并給出徑向 分 布函 數g(r)=g(r,Reλ,St)， 這 里的g(r)不 僅 包含Stokes數，還包含了以泰勒微尺度為特征尺度的雷諾數（U′為湍流強度）。當St=1時，增強因子最大，當St→0 或St→∞時，增強因子趨近于1。此外，雷諾數Reλ越大，增強因子越大。

Zhou等[71]得到雙分散顆粒情況下任意兩種尺度顆粒的凝并核函數：

Saffman和Turner[56]提出適用于剪切條件下的湍流凝并核，Chun等[74]基于該凝并核，得到在St→0情況下具有顆粒聚集效應的湍流凝并核：

式中c1∝St2， γ≈1， 在St→0 且 η/dij?1的情況下，顆粒聚集的增強效應接近于c0（一個待定的匹配系數）。

上述研究都是基于顆粒的拉格朗日描述用DNS來求解流場，這在慣性顆粒系統中有效[75-76]，但對零慣性顆粒（忽略顆粒質量）而言則無效。從物理機制看，零慣性顆粒有良好的跟隨性，顆粒與流體的速度幾乎一致，湍流場的渦運動無法因離心效應而導致顆粒的聚集。因此，需要采用其他方法描述零慣性顆粒的濃度脈動。

3.2 湍流脈動對凝并生長效應的影響

顆粒的生長除了凝積外，還可以通過凝并。凝并使顆粒粒徑增大、濃度減少。定義粒徑增長率為：

式中dg和vg是顆粒的幾何平均直徑和幾何平均體積。粒徑增長率可用來綜合衡量顆粒成核、凝積生長和凝并導致顆粒粒徑增長的效果，在單獨考慮凝并作用時，也可以衡量凝并對顆粒粒徑的增長作用。

在直接數值模擬方法中，對包含完整脈動信息的顆粒數密度場n的模擬，可以得到粒徑的瞬時增長率Ω。對采用大尺度濾波的數密度場 〈n〉模擬，得到的時濾波增長率 ΩF，則亞網格尺度（SGS）作用增長率為Garrick[36]采用DNS和LES方法研究混合層湍流，通過矩方法使包含顆粒凝并源項的方程得以封閉后求解，得到了瞬時增長率 Ω和SGS尺度作用增長率 ΩSGS，并建立SGS尺度作用增長率與流場剪切旋轉效應和湍流場的標量耗散率之間的關系（如圖4所示）。

圖4 SGS尺度粒徑增長效應Ω SGS 的概率密度函數分布[36]Fig. 4 PDF of SGS particle size growth Ω SGS[36]

在湍流場中，被動標量M1的耗散率定義為：

式中M1是n(v)的內部一階矩，與顆粒的總質量濃度相關，不受顆粒凝并作用的影響。由圖4（b）可知，在質量濃度耗散率為零和最大值處， ΩSGS的PDF分布比較對稱，且集中分布在Q=0附近的區域。在零和最大值處之間， ΩSGS的PDF分布非常分散且偏向負值。可見，由湍流脈動導致的顆粒濃度脈動抑制了顆粒的凝并，這與Rigopoulos[25]的研究結論相符。綜合顆粒的成核、凝積生長和凝并這幾個因素，湍流脈動會抑制顆粒的增長[14,18,37]。

3.3 顆粒凝并項與破碎項增強因子的重新定義

于是，可用gnn來對顆粒凝并項進行描述，即：

該式中還需確定gnn。Yang等[77]認為對零慣性顆粒而言，gnn可表示為流場湍動能與總動能之比：

式中k是湍動能，他們基于該式對圓管湍流場進行了數值模擬，得到的結果與實驗結果吻合。

定義新的增強因子：

令g3=gnn+2gΓn+gΓnn，則再代入凝并項后得[78]：

則新的增強因子定義為gq

此處沿用并定義gnn、g3和gq等增強因子出于兩方面考慮，一是基于顆粒聚集效應的研究[73-74,79]，在極低顆粒慣性情況下，顆粒凝并的增強受顆粒多分散性的作用很小，則二是若忽略顆粒多分散性的作用，則gnn近似于顆粒數密度脈動的均方根，這便于對湍流脈動凝并項方程的封閉求解，甚至可采用Elghobashi[64]提出的k-ε-g方程求解。關于g3增強因子中出現的gΓn、gΓnn，已有前人進行過研究[35]。

4 結束語

微納米顆粒兩相湍流存在于自然現象和大量的實際應用中，研究其規律具有重要的實際意義。微納米顆粒兩相流與大顆粒兩相流相比具有其特殊性和復雜性，對其研究具有重要的科學意義。

在微納米顆粒兩相湍流中，湍流導致的流場速度脈動、顆粒濃度脈動、溫度脈動以及飽和度脈動會對顆粒相的成核、凝積生長、碰撞凝并、破碎產生影響。本文在給出顆粒一般動力學方程的基礎上，分析了方程中因湍流脈動導致的與顆粒成核、生長、凝并、破碎相關的脈動量關聯項，介紹了對這些關聯項進行模化從而使方程得以封閉的研究途徑和相關的研究結果。

湍流脈動對顆粒成核、凝積生長、碰撞凝并、破碎的影響已受到人們的重視并開展了研究，但仍有一些方面的問題有待將來進一步探討。

1）顆粒聚集使顆粒碰撞的幾率增加，從而導致顆粒凝并的增強。而在第3.3節的研究結論中，湍流脈動會抑制顆粒的凝并，這與顆粒聚集導致顆粒凝并增強的結論似乎存在矛盾，因為湍流脈動使得顆粒數密度脈動，而數密度脈動的結果也有可能導致顆粒聚集。所以有必要對湍流脈動對顆粒凝并的影響進行深入的研究。

2）微納顆粒兩相湍流場的實驗研究具有一定的難度。湍流場尤其是充分發展的湍流場具有空間和時間多尺度的特征。受多尺度的影響，在測量顆粒的取樣頻率還很低的情況下，難以全面獲得微納尺度顆粒在流場中的運動信息[80]，有必要在實驗儀器和實驗技術方面做進一步的研究。

3）在納顆粒兩相湍流場中，對顆粒一般動力學方程進行雷諾平均后會出現與顆粒成核、生長、凝并、破碎相關的脈動量關聯項，這些項使得平均后的方程不封閉，要使方程封閉，就要對這些項進行模化。目前，由模化后得到的數值模擬結果與實驗結果如顆粒燃燒[81-83]、成核[13]、凝并[84]、破碎彌散[85-86]還存在一定的差異，原因是數值模擬中可以有針對性地對成核、生長、凝并、破碎中的某一項進行，而實驗的結果往往是這幾項中若干項的綜合結果。所以有必要進一步對顆粒一般動力學方程中未知項的模化開展有針對性的實驗研究。