湍流邊界層中顆粒物擴(kuò)散的雙拉格朗日模擬計算

卓　俊，黃流興，牛勝利，謝紅剛

(西北核技術(shù)研究所,西安710024)

?

湍流邊界層中顆粒物擴(kuò)散的雙拉格朗日模擬計算

卓俊，黃流興，牛勝利，謝紅剛

(西北核技術(shù)研究所,西安710024)

采用拉格朗日方法模擬計算了湍流邊界層中氣相粒子和顆粒物的運(yùn)動。其中，氣相粒子在流場中按平均風(fēng)速輸送，用一系列隨機(jī)位移，模擬湍流對擴(kuò)散的影響，分析了顆粒物在流場中的受力情況，建立了求解顆粒物擴(kuò)散軌跡的牛頓運(yùn)動方程。用該方法模擬平板邊界層內(nèi)顆粒物的擴(kuò)散行為，給出了顆粒物受力及在邊界層內(nèi)的運(yùn)動軌跡，并將計算的質(zhì)量通量與文獻(xiàn)中的實測值進(jìn)行了比較。結(jié)果表明：計算值與實測值基本吻合，初步驗證了本文方法在模擬顆粒物擴(kuò)散時的可靠性。

氣固兩相流；顆粒物；湍流；擴(kuò)散；雙拉格朗日

顆粒流體系統(tǒng)在自然界中普遍存在，與生態(tài)環(huán)境、工業(yè)過程及人們?nèi)粘Ｉ畹母鱾€方面密切相關(guān)。顆粒流體系統(tǒng)涉及顆粒物輸運(yùn)以及流場與周圍環(huán)境的相互作用，準(zhǔn)確預(yù)測顆粒物擴(kuò)散是環(huán)境工程和工農(nóng)業(yè)生產(chǎn)中的重要課題之一[1-3]。同時，湍流流動是自然界常見的現(xiàn)象，在多數(shù)工程問題中流體的流動往往處于湍流狀態(tài)，湍流流動特征在工程中占有重要的地位，湍流研究一直受到研究者的高度重視。

目前，預(yù)測粒子在湍流場中擴(kuò)散的數(shù)值計算方法主要有歐拉-歐拉方法和歐拉-拉格朗日方法兩種[4-5]。歐拉-歐拉方法，也稱為顆粒相擬流體模型，已廣泛應(yīng)用于稠密顆粒相的氣固流動體系中。該模型把顆粒相處理為具有連續(xù)介質(zhì)特性的、與氣相相互滲透的擬流體，氣相和顆粒相都在歐拉坐標(biāo)下進(jìn)行求解，仿照單相流動對顆粒湍流脈動進(jìn)行模擬，應(yīng)用顆粒動力學(xué)理論及分子運(yùn)動論使方程組封閉。顆粒相擬流體模型最大的優(yōu)點是可以全面地考慮顆粒相的輸運(yùn)特性，并進(jìn)行大規(guī)模的工程問題計算。但其中的顆粒連續(xù)性假設(shè)，削弱了氣固兩相流中非均勻結(jié)構(gòu)的真實性，計算中無法體現(xiàn)顆粒的大小、密度等物理性質(zhì)的區(qū)別，模擬結(jié)果也無法揭示顆粒復(fù)雜的變化經(jīng)歷[4]。歐拉-拉格朗日方法將流體相處理為連續(xù)相，在歐拉坐標(biāo)系下，建立維納-斯托克斯方程組求解連續(xù)相的流動、傳熱及反應(yīng)特性，而在拉格朗日坐標(biāo)系下，應(yīng)用牛頓第二定律跟蹤求解流場中的每一個離散粒子的運(yùn)動軌跡，以反映整個離散顆粒場；連續(xù)相-離散相的相互作用服從牛頓第三定律，以源相添加于各自的求解方程中實現(xiàn)相間的作用耦合。在該模型中，如何合理描述氣相湍流場對顆粒運(yùn)動的影響是一個主要難題。雖然參照氣相對湍流擴(kuò)散的拉格朗日描述，引入了顆粒湍流運(yùn)動的拉格朗日時間尺度概念，但不同尺度顆粒物的拉格朗日時間尺度計算方法及其計算結(jié)果的可靠性，都還有待進(jìn)一步研究。

為克服歐拉-拉格朗日方法在模擬計算中遇到的困難，近年，有學(xué)者提出了氣相和固相均用拉格朗日方法描述的雙拉格朗日方法，用來模擬計算顆粒相在湍流中的擴(kuò)散過程[6-7]。該方法適用于顆粒濃度很低，可以忽略氣相中由顆粒引起的擾流情況，即稀相流情形。根據(jù)氣固兩相系統(tǒng)中氣相運(yùn)動的特點，將兩相流動中顆粒的運(yùn)動過程分解為受沖力支配的瞬時碰撞運(yùn)動及受流體曳力控制的懸浮運(yùn)動，從而建立了顆粒運(yùn)動模型。該模型中，流體和顆粒的運(yùn)動規(guī)律均用拉格朗日方法求解。

本文介紹了雙拉格朗日方法的物理模型，建立了基于雙拉格朗日方法的氣固兩相流數(shù)值計算模型，用編制的計算機(jī)程序，模擬計算了無限平板邊界層中固體顆粒的擴(kuò)散過程，并將計算結(jié)果與文獻(xiàn)中的實驗數(shù)據(jù)進(jìn)行了對比，從而對本文方法的正確性進(jìn)行了初步驗證。

1物理模型

雙拉格朗日方法中，用大量標(biāo)記粒子表征氣相和固相粒子的擴(kuò)散，如圖1所示。假設(shè)在某時刻t，氣相粒子和顆粒物在流場中處于相同的位置，即Xg(t)=Xp(t)，記該點為O，并假設(shè)此時氣相粒子和顆粒物的瞬時速度已知，分別為Vg(t)和Vp(t)，在經(jīng)過一段時間間隔Δt后，氣相粒子和顆粒物分別運(yùn)動至F點和P點，其坐標(biāo)分別為

(1)

(2)

在t+Δt時刻，粒子的瞬時速度Vp(t+Δt)可由顆粒運(yùn)動方程計算得到：

(3)

其中，mp為所跟蹤固體顆粒的質(zhì)量；Fi表示固體顆粒在流場中所受的作用力。在計算固體顆粒在P點附近的受力大小時，需要預(yù)先知道同在P點處的氣相瞬時速度VP,g(t+Δt)。在已知O點氣相粒子運(yùn)動速度VO,g(t)的情況下，可通過氣相湍流速度的拉格朗日時間自相關(guān)關(guān)系給出氣相粒子在F點的瞬時速度VF,g(t+Δt)，再通過氣相瞬時運(yùn)動速度的歐拉空間關(guān)聯(lián)，求得P點的氣相速度。

圖1雙拉格朗日方法示意圖Fig.1Sketch of Lagrangian-Lagrangian simulation approach

2計算方法

在兩相流模擬中，氣相粒子在t時刻的位置為

(4)

其中，vg,x(x,y,z,t-Δt)、vg,y(x,y,z,t-Δt)和vg,z(x,y,z,t-Δt)分別為t-Δt時刻氣相粒子在x、y、z坐標(biāo)方向的速度分量。在湍流場中，氣相的速度由平均速度和脈動速度組成，即

(5)

(6)

其中，σx、σy和σz分別為脈動風(fēng)速在x、y、z坐標(biāo)方向的標(biāo)準(zhǔn)差；ξ是一個由計算機(jī)產(chǎn)生的標(biāo)準(zhǔn)正態(tài)分布隨機(jī)數(shù)；RL,x、RL,y和RL,z分別為x、y、z坐標(biāo)方向的拉格朗日自相關(guān)系數(shù)：

(7)

其中，TL,x、TL,y和TL,z分別是x、y、z坐標(biāo)方向的拉格朗日相關(guān)時間尺度。

另外，在同一時刻不同位置上的湍流速度，可以采用與式(6)類似的歐拉關(guān)聯(lián)計算：

(8)

其中，RE，x、RE，y和RE，z分別為x、y、z坐標(biāo)方向的歐拉相關(guān)系數(shù)，由歐拉時間尺度通過與式(7)類似的關(guān)系計算得到。拉格朗日相關(guān)時間尺度和歐拉相關(guān)時間尺度之比β可表示為

(9)

其中，i為湍流強(qiáng)度，可表示為

(10)

顆粒在氣流中運(yùn)動時，主要受到重力G、黏性阻力FD、Basset力FB、Magnus力FM、附加質(zhì)量力Fmv、Saffman升力FS和壓力梯度力Fp的作用。假設(shè)顆粒為不變形、無相互碰撞的球形顆粒，其存在對流體相無影響，那么，顆粒在邊界層中運(yùn)動的受力關(guān)系為

(12)

基于上述方法，編制了顆粒物在氣流中運(yùn)動的氣固兩相流模擬程序。該程序能夠計算不同粒徑顆粒物在已知平均流場時空分布的氣流中的運(yùn)動情況，程序流程圖如圖2所示。其中，Vg(Xg,t+Δt)和Vg(Xp,t+Δt)分別為氣流在Xg(t+Δt)和Xp(t+Δt)處的運(yùn)動速度。

圖2雙拉格朗日程序流程圖Fig.2The flow chart of Lagrangian-Lagrangian simulation code

3結(jié)果及分析

為了驗證上述方法的正確性，對1997年Tanière等所做的顆粒物擴(kuò)散實驗進(jìn)行了模擬計算。Tanière實驗研究了玻璃粒子和聚氯乙烯(PVC)粒子在平板邊界層中的擴(kuò)散過程[10]。該實驗中氣相自由速度Ve=10.6 m·s-1, 摩擦速度V*= 0.4 m·s-1，玻璃粒子直徑dglass=60 μm，密度ρglass=2 500 kg·m-3，PVC粒子直徑dPVC=130 μm，

密度ρPVC=1 430 kg·m-3，顆粒物在離平板端點3.17 m處從平板底部進(jìn)入流場，釋放初始速度vy0為2.0V*～2.5V*，入射角度α為15°～45°。根據(jù)實驗測量參數(shù)，實驗中的一些物理量分布如下[11]：

邊界層厚度：

(13)

邊界層速度分布：

(14)

拉格朗日相關(guān)時間尺度：

(15)

邊界層內(nèi)湍流速度的標(biāo)準(zhǔn)差：

(16)

此外，計算中涉及的空氣密度ρa(bǔ)為1.19 kg·m-3，空氣動力黏度μ為1.81×10-5Pa·s，重力加速度g為9.81 m·s-2。計算中粒子碰壁過程按彈性碰撞處理。

表1為顆粒物在離平板端點5 m處的受力情況。計算Magnus力FM時，顆粒的轉(zhuǎn)速取1 000 r·s-1。由于Saffman升力FS和壓力梯度力Fp的大小與顆粒到邊界的距離密切相關(guān)，因此表中同時給出了粒子所在高度為0.01δ(x)，0.1δ(x)和0.5δ(x)時FS和Fp的值。從表1可以看出，在上述給定的物理特征參數(shù)下，玻璃和PVC顆粒在平板邊界層中的運(yùn)動行為，主要受氣流黏性阻力和Basset力的支配，但由于重力方向始終豎直向下，重力累積效應(yīng)導(dǎo)致的粒子沉降也不可忽略。另外，當(dāng)顆粒比較接近邊界位置時，氣流在垂直方向有較大的速度梯度，此時，顆粒受到的Saffman升力FS和壓力梯度力Fp會明顯增大，且這兩種力的方向始終豎直向上，F(xiàn)S約比Fp大1個數(shù)量級，因此，在邊界層中，Saffman力是使顆粒物不發(fā)生沉降的主要作用力。

表1顆粒物受到各種力的大小

圖3為玻璃粒子和PVC粒子的運(yùn)動軌跡。從圖中可以看出，兩種粒子在邊界層底部做跳躍運(yùn)動，只有個別粒子能夠進(jìn)入主流區(qū)。PVC粒子的運(yùn)動軌跡類似在重力作用下的拋物線，而玻璃粒子由于質(zhì)量較小，慣性效應(yīng)較弱，運(yùn)動更容易受湍流波動的影響，因此軌跡明顯偏離只有重力作用時的拋物線運(yùn)動，表現(xiàn)出較強(qiáng)的隨機(jī)性。

(a)Glass

(b)PVC圖3粒子軌跡Fig.3Particle traces

跟蹤模擬計算了10 000個粒子在邊界層中的運(yùn)動軌跡，圖4給出了離平板端點5.06 m處顆粒物質(zhì)量通量Φm隨離開壁面法向無量綱距離y+的變化情況，并將模擬計算結(jié)果與文獻(xiàn)[10]的實測結(jié)果進(jìn)行了比較。其中，y+可表示為

(17)

從圖4可以看出，計算結(jié)果與實驗結(jié)果符合較好，變化趨勢也基本一致，初步驗證了本文方法在顆粒物擴(kuò)散計算中的可行性。但圖4中的計算結(jié)果與實驗值還存在一定的偏差。引起偏差的原因可能是在本文的計算中忽略了顆粒間的碰撞。此外，由于目前對湍流的認(rèn)識不足，計算中采用的湍流特征參數(shù)與實驗中的實際情況可能存在偏差。

(a)Glass

(b)PVC圖4質(zhì)量通量隨y方向無量綱距離y+的變化Fig.4Mass flux vs. y+

4小結(jié)

建立了一種對氣相和固體顆粒均采用拉格朗日方法模擬的雙拉格朗日粒子擴(kuò)散計算模型。在計算氣相流動時，用拉格朗日自相關(guān)關(guān)系計算相同粒子在不同時刻的湍流速度，用歐拉關(guān)聯(lián)描述不同位置上的氣相粒子在同一時刻的湍流速度關(guān)系；使用牛頓運(yùn)動方程描述固體顆粒物在流體中的運(yùn)動行為，考慮了顆粒物在流體中受到的黏性阻力、Saffman升力、Basset力及重力的影響。依據(jù)該模型編寫了計算程序，模擬了顆粒物在平板邊界層的運(yùn)動，并將計算結(jié)果與文獻(xiàn)的實驗結(jié)果進(jìn)行了比較。分析認(rèn)為：計算結(jié)果與文獻(xiàn)實測結(jié)果吻合良好。

[1]李靜海, 歐陽潔, 高士秋, 等. 顆粒流體復(fù)雜系統(tǒng)的多尺度模擬[M]. 北京: 科學(xué)出版社, 2005.(LI Jing-hai, OUYANG Jie, GAO Shi-qiu, et al. Multi-Scale Simulation of Particle-Fluid Complex Systems[M]. Beijing: Science Press, 2005.)

[2]WANG X M, DONG Z B, ZHANG J W, et al. Relations between morphology, air flow, sand flux and particle size on transverse dunes, Taklimakan sand sea, China[J]. Earth Surface Processes and Landforms: The Journal of the British Geomorphological Research Group, 2002, 27(5): 515-526.

[3]危衛(wèi), 魯錄義, 顧兆林. 風(fēng)沙運(yùn)動的電場-流場耦合模型及氣固兩相流數(shù)值模擬[J]. 物理學(xué)報, 2012，61(15): 158301. (WEI Wei, LU Lu-yi, GU Zhao-lin. Modeling and simulation of electrification of wind-blown-sand two-phase flow[J]. Acta Phys Sin, 2012, 61(15): 158301.)

[4]馬銀亮. 高濃度氣固兩相流的數(shù)值模擬研究[D]. 杭州: 浙江大學(xué), 2001.(MA Yin-liang. Study on numerical simulation of dense gas-particle two-phase flow[D]. Hangzhou: Zhejiang University, 2001.)

[5]劉向軍, 石磊, 徐旭常. 稠密氣固兩相流歐拉-拉格朗日的研究現(xiàn)狀[J]. 計算力學(xué)學(xué)報, 2007, 24(2): 166-172.(LIU Xiang-jun, SHI Lei, XU Xu-chang. Activities of dense particle-gas two-phase flow modeling in Eulerian-Lagrangian approach[J]. Chinese Journal of Computational Mechanics, 2007, 24(2): 166-172.)

[6]昌澤舟, BERLEMONT A, GOUESBET G. 氣固兩相流中顆粒彌散的拉格朗日模擬[J]. 計算力學(xué)學(xué)報, 2001, 18(3): 284-289.(CHANG Ze-zhou, BERLEMONT A, GOUESBET G. Lagrangian simulation of particle dispersion in gas-solid turbulent flows[J]. Chinese Journal of Computational Mechanics, 2001, 18(3): 284-289.)

[7]CHEN B, WANG C, WANG Z W, et al. Investigation of gas-solid flow across circular cylinders with discrete vortex method[J]. Applied Thermal Engineering, 2009, 29(8/9): 1 457-1 466.

[8]THOMSON D J. Criteria for the selection of stochastic models of particle trajectories in turbulent flows[J]. J Fluid Mech, 1987, 180: 520-556.

[9]CLIFT R, GRACE J R, WEBER M E. Bubbles, Drops and Particles[M]. New York: Academic Press, 1978.

[11]郭宇, 崔桂香, 許春曉, 等. 湍流邊界層中重粒子彌散的隨機(jī)模型[J]. 計算物理, 2004, 21(6): 515-522.(GUO Yu, CUI Gui-xiang, XU Chun-xiao, et al. A stochastic model of heavy particle dispersion in turbulent boundary layer[J]. Chinese Journal of Computational Physics, 2004, 21(6): 515-522.)

Lagrangian-Lagrangian Simulation of Gas-Solid Two-Phase Flow in Turbulent Boundary Layer

ZHUO Jun,HUANG Liu-xing,NIU Sheng-li,XIE Hong-gang

(Northwest Institute of Nuclear Technology,Xi’an710024,China)

A Lagrangian-Lagrangian simulation approach of gas-solid two-phase flow in boundary layer with turbulence is presented, and the reliability verification of the approach is also provided in this paper. In the Lagrangian-Lagrangian simulation approach, both gas and solids are simulated in Lagrangian model. The gas is transported under mean wind, meanwhile, the influence of turbulence is dealt with a series of random movements of gas. The velocities of turbulence sequence in time are simulated by Lagrangian relationship, and velocities of turbulence at different positions are determined by Eular relationship. The particle traces are determined by solving Newtonian equations of motion. The forces such as drag forces, Magnus forces, Basset forces and gravity are considered when particles move in gas flow. The proposed model is applied to simulating a horizontal turbulent boundary layer of air carrying heavy solid particles. The values of forces imposed on particles and the images of particle traces during transportation are obtained. The calculated results of mass flux of particles have been compared with the experimental data, and they are in good agreement.

gas-solid two-phase flow;solid particles;turbulence;diffusion;Lagrangian-Lagrangian

2016-07-04；

2016-07-18

卓俊(1982-)，男，浙江龍泉人，助理研究員，碩士，主要從事放射性核素大氣擴(kuò)散數(shù)值模擬研究。

E-mail:zhuojun@nint.ac.cn

O359

A

2095-6223(2016)031004(6)