郭孝武, 凡鳳仙, 胡曉紅, 蘇明旭
(1.上海理工大學(xué) 能源與動(dòng)力工程學(xué)院,上海 200093;2.上海理工大學(xué) 上海市動(dòng)力工程多相流動(dòng)與傳熱重點(diǎn)實(shí)驗(yàn)室,上海 200093)
?
特征參數(shù)對(duì)DSMC方法模擬聲場(chǎng)中顆粒碰撞的影響
郭孝武1,2, 凡鳳仙1,2, 胡曉紅1,2, 蘇明旭1,2
(1.上海理工大學(xué) 能源與動(dòng)力工程學(xué)院,上海 200093;2.上海理工大學(xué) 上海市動(dòng)力工程多相流動(dòng)與傳熱重點(diǎn)實(shí)驗(yàn)室,上海 200093)
直接從顆粒受力和運(yùn)動(dòng)出發(fā),基于直接模擬蒙特卡洛(DSMC)方法建立顆粒碰撞模型,模擬聲場(chǎng)中顆粒的碰撞過(guò)程,通過(guò)改變模擬條件,探討DSMC方法中的特征參數(shù)(取樣顆粒的數(shù)目權(quán)重、網(wǎng)格數(shù)目、時(shí)間步長(zhǎng))對(duì)顆粒碰撞率和計(jì)算時(shí)間的影響.結(jié)果表明:頻率越大,顆粒容積份額變化越迅速,顆??臻g分布越不均勻;數(shù)目權(quán)重的增加對(duì)顆粒碰撞率影響較小,而對(duì)計(jì)算時(shí)間的影響顯著,使其迅速減少;網(wǎng)格數(shù)目增加,碰撞率降低,計(jì)算時(shí)間則先迅速降低,隨后在低頻時(shí)基本不變,高頻時(shí)略有上升.研究還發(fā)現(xiàn):隨著碰撞時(shí)間步長(zhǎng)的增加,低頻聲場(chǎng)中碰撞率單調(diào)增加,高頻聲場(chǎng)中碰撞率先增加而后出現(xiàn)波動(dòng);碰撞時(shí)間步長(zhǎng)的增加將引起計(jì)算時(shí)間減少,減少量在碰撞步長(zhǎng)較小時(shí)最為明顯.
聲場(chǎng); 顆粒; 碰撞; 直接模擬蒙特卡洛方法; 特征參數(shù)
目前,用于聲凝并的數(shù)值模擬方法可歸納為區(qū)域算法、矩量法和直接模擬蒙特卡洛(direct simulation Monte Carlo,DSMC)方法[8].其中,DSMC方法立足于氣體分子運(yùn)動(dòng)論,將所計(jì)算的實(shí)際顆粒場(chǎng)用取樣顆粒場(chǎng)進(jìn)行置換,對(duì)每一個(gè)取樣顆粒賦予一個(gè)數(shù)目權(quán)重(即取樣顆粒所代表的真實(shí)顆粒數(shù)目),跟蹤取樣顆粒的運(yùn)動(dòng)軌跡,通過(guò)概率的方法判斷碰撞是否發(fā)生,并對(duì)碰撞后的顆粒凝并事件進(jìn)行處理,該方法可以方便地分析多分散顆粒凝并的動(dòng)態(tài)過(guò)程以及粒徑演變規(guī)律[9-12].2006年,Sheng等[9]利用常體積Monte Carlo方法在初始顆粒及團(tuán)聚體均為球形的前提下,研究了行波聲場(chǎng)中液態(tài)PM2.5的聲凝并,發(fā)現(xiàn)數(shù)值模擬結(jié)果在一些情況下與實(shí)驗(yàn)吻合較好,而當(dāng)顆粒初始粒徑分布發(fā)生改變時(shí),則與實(shí)驗(yàn)存在明顯差異;2007年,他們引入分形維數(shù)以描述顆粒團(tuán)聚體的形狀和結(jié)構(gòu),將上述聲凝并模型推廣到固體顆粒的聲凝并中[10].與Sheng等[9-10]在顆粒通用動(dòng)力學(xué)方程和凝并核函數(shù)的基礎(chǔ)上建立聲凝并模型不同,凡鳳仙等[11-13]從顆粒運(yùn)動(dòng)方程出發(fā),利用DSMC方法研究PM2.5在行波聲場(chǎng)中的碰撞特性[11]、在駐波聲場(chǎng)中的碰撞和凝并規(guī)律[12-13].顆粒碰撞是其發(fā)生凝并的基礎(chǔ)和前提,然而,在利用DSMC方法處理顆粒碰撞時(shí),3個(gè)特征參數(shù)(時(shí)間步長(zhǎng)、取樣顆粒的數(shù)目權(quán)重、計(jì)算區(qū)域內(nèi)的網(wǎng)格數(shù)目)需要人為設(shè)定,這些特征參數(shù)的取值將影響到計(jì)算精度[14-16].時(shí)間步長(zhǎng)與取樣顆粒的數(shù)目權(quán)重越小,計(jì)算越精確,然而隨之出現(xiàn)的是計(jì)算代價(jià)的增加;若在一個(gè)網(wǎng)格范圍內(nèi)判斷顆粒碰撞的發(fā)生,則計(jì)算區(qū)域內(nèi)的網(wǎng)格數(shù)目也將對(duì)計(jì)算結(jié)果和計(jì)算代價(jià)帶來(lái)影響.因此,本文基于DSMC方法,在不同特征參數(shù)條件下對(duì)聲場(chǎng)中顆粒的碰撞過(guò)程進(jìn)行數(shù)值模擬,確定特征參數(shù)對(duì)顆粒碰撞率和計(jì)算時(shí)間的影響,為利用DSMC方法準(zhǔn)確、高效地研究PM10,PM2.5的聲凝并提供參考.
為利用DSMC方法研究聲場(chǎng)中顆粒的碰撞過(guò)程,作出如下簡(jiǎn)化假設(shè):
a. 聲場(chǎng)為一維平面駐波聲場(chǎng),聲波波動(dòng)方向?yàn)樗椒较?
b. 由于煙氣的物性參數(shù)和空氣類(lèi)似,且其壓力不太高,將氣體介質(zhì)視為理想空氣;
c. 顆粒為剛性球體,忽略顆粒的轉(zhuǎn)動(dòng);
d. 由于聲輻射壓力的作用效果極其微弱[17],不考慮聲輻射壓力;
e. 顆粒碰撞在各網(wǎng)格內(nèi)進(jìn)行,且為二元碰撞,為著重考慮顆粒的碰撞,暫不考慮其凝并,認(rèn)為顆粒碰撞后彈開(kāi);
f. 在顆粒碰撞過(guò)程中,顆粒發(fā)生滑移時(shí)所受摩擦力遵守庫(kù)侖摩擦定律.
1.1 聲波波動(dòng)方程
由Navier-Stokes方程,可推導(dǎo)出無(wú)旋、無(wú)粘流體中x向駐波聲場(chǎng)的波動(dòng)方程為
(1)
式中:ufx為聲波引起的流體介質(zhì)振動(dòng)速度;x為位置坐標(biāo);t為時(shí)間;ua為速度振幅;k為波數(shù),k=ω/c,c為聲速;ω=2πf,f為聲場(chǎng)頻率.
習(xí)慣上常用聲壓級(jí)和頻率描述聲場(chǎng),聲壓級(jí)與速度振幅ua的關(guān)系為
(2)
式中:L為聲壓級(jí),dB;Pr為參考聲壓,Pr=210-5Pa.
1.2 顆粒運(yùn)動(dòng)方程
設(shè)重力方向?yàn)閦向,與x和z垂直的方向?yàn)閥向,流體攜帶顆粒沿y向以一定速度通過(guò)聲場(chǎng)空間.根據(jù)牛頓第二定律,顆粒運(yùn)動(dòng)方程可寫(xiě)為
其中:式(3a)等號(hào)右邊4項(xiàng)依次為Stokes力的x向分量、壓力梯度力、虛擬質(zhì)量力、Basset力;式(3b)等號(hào)右邊為Stokes力的y向分量;式(3c)等號(hào)右邊3項(xiàng)分別為Stokes力的z向分量、重力、浮力.式中:mp為顆粒質(zhì)量;mf為與顆粒等體積的流體質(zhì)量;ufy和ufz分別為流體速度的y,z向分量;upy與upz分別為顆粒速度的y,z向分量;ρp為顆粒密度;g為重力加速度;t′為時(shí)間變量;Cc為Cunningham修正系數(shù),其表達(dá)式為[13]
式中,Kn為Knudsen數(shù),Kn=2lm/dp,lm為氣體分子平均自由程,m.
1.3 描述顆粒碰撞的DSMC方法
對(duì)含有大量PM10,PM2.5的氣固懸浮體系施加駐波聲場(chǎng),顆粒間的碰撞不可避免.為了準(zhǔn)確地描述聲場(chǎng)中顆粒的行為規(guī)律,理想的方法是跟蹤到每一個(gè)顆粒,通過(guò)顆粒運(yùn)動(dòng)軌道來(lái)判斷顆粒碰撞,但當(dāng)顆粒量很大時(shí),這種方法的計(jì)算量將是非常驚人的.由于顆粒數(shù)目龐大,利用顆粒軌道判斷碰撞的方法難以勝任.解決這一問(wèn)題的一個(gè)行之有效的途徑是采用DSMC方法.由于DSMC方法處理的是遠(yuǎn)低于真實(shí)顆粒數(shù)目的取樣顆粒,從而大大降低了計(jì)算量,特別適于在較大規(guī)模顆粒量和較大區(qū)域范圍內(nèi)對(duì)小粒徑顆粒、較稀的氣固兩相流動(dòng)進(jìn)行模擬研究.在DSMC方法的前提下,不同研究者的處理和計(jì)算方法不盡相同,較常用的是將計(jì)算區(qū)域劃分成若干個(gè)網(wǎng)格,在同一網(wǎng)格內(nèi)判斷顆粒之間是否發(fā)生碰撞,本文也采用該方法.
(5)
取樣顆粒i與同一網(wǎng)格中其他顆粒的總碰撞概率Pi可表示為
式中,N為取樣顆粒i所在網(wǎng)格的取樣顆??倲?shù).
采用修正的Nanbu方法[19]判斷顆粒碰撞的發(fā)生,即對(duì)于任意取樣顆粒i,在滿(mǎn)足Pi<1的前提下,利用[0,1)區(qū)間上均勻分布的隨機(jī)數(shù)R,選擇同一網(wǎng)格內(nèi)的候選被碰取樣顆粒j為
j=int[RN]+1
(7)
式中,int[RN]表示RN的整數(shù)部分.如果R滿(mǎn)足
(8)
則認(rèn)為取樣顆粒i與j發(fā)生碰撞.此時(shí),保持顆粒位置不變,顆粒碰撞后速度根據(jù)動(dòng)量定理計(jì)算[19],即
(9)
(10)
(11)
(12)
(13)
(14)
(15)
2.1 邊界條件與計(jì)算參數(shù)
取λ2k×1 mm×1 mm的三維空間為計(jì)算區(qū)域,其中λ2k表示頻率f=2 kHz時(shí)對(duì)應(yīng)的聲波波長(zhǎng),沿x向?qū)⒂?jì)算區(qū)域劃分為若干個(gè)大小相等的網(wǎng)格.由聲波的周期性和本文采用的聲波頻率(f=2,4,6 kHz)可知,對(duì)x,y,z這3個(gè)方向的邊界均可采用周期性邊界條件[13]進(jìn)行處理,對(duì)各個(gè)取樣顆粒采用相同的數(shù)目權(quán)重,主要計(jì)算參數(shù)見(jiàn)表1.其中:T為氣體溫度;p為氣體靜態(tài)壓力;np為顆粒數(shù)目濃度;Δt為求解顆粒運(yùn)動(dòng)的時(shí)間步長(zhǎng).本文通過(guò)試算確定Δt,表1中給出的Δt滿(mǎn)足3種計(jì)算頻率下顆粒運(yùn)動(dòng)的計(jì)算精度要求.
表1 數(shù)值計(jì)算參數(shù)
2.2 顆粒運(yùn)動(dòng)方程的求解
采用定步長(zhǎng)四階Runge-Kutta算法對(duì)式(3)進(jìn)行求解.其中,Basset力與顆粒經(jīng)歷的運(yùn)動(dòng)過(guò)程有關(guān).將流體速度和顆粒速度基于經(jīng)歷的時(shí)間步離散,可實(shí)現(xiàn)Basset力的數(shù)值求解.例如,在時(shí)間區(qū)間0~t3離散為0~t1,t1~t2,t2~t3的情況下,如果這3個(gè)時(shí)間區(qū)間對(duì)應(yīng)的流體與顆粒速度增量分別為Δufx0,Δufx1,Δufx2與Δupx0,Δupx1,Δupx2,則t3時(shí)刻顆粒受到的Basset力FB可按下式求解[20]:
(16)
顆粒運(yùn)動(dòng)軌跡可由各時(shí)間步長(zhǎng)內(nèi)顆粒的位移疊加得到,根據(jù)二階隱式Adams插值算法,顆粒的位移Sp可表示為
式中,up為顆粒速度.
2.3 數(shù)值計(jì)算流程
圖1 計(jì)算流程
3.1 聲場(chǎng)中顆粒容積份額的演變
圖2 不同頻率的聲場(chǎng)中顆粒容積份額的演變
3.2 特征參數(shù)對(duì)碰撞率和計(jì)算時(shí)間的影響
3.2.1 數(shù)目權(quán)重對(duì)碰撞率和計(jì)算時(shí)間的影響
3.2.2 網(wǎng)格數(shù)目對(duì)碰撞率和計(jì)算時(shí)間的影響
圖3 數(shù)目權(quán)重對(duì)碰撞率和計(jì)算時(shí)間的影響
圖4 網(wǎng)格數(shù)目對(duì)碰撞率和計(jì)算時(shí)間的影響
3.2.3 碰撞時(shí)間步長(zhǎng)對(duì)碰撞率和計(jì)算時(shí)間的影響
圖5 碰撞時(shí)間步長(zhǎng)對(duì)碰撞率和計(jì)算時(shí)間的影響
綜合考慮顆粒受到的Stokes力、壓力梯度力、虛擬質(zhì)量力、Basset力、重力和浮力,建立顆粒運(yùn)動(dòng)模型,基于DSMC方法模擬聲場(chǎng)中顆粒的碰撞過(guò)程,獲得顆粒碰撞率和計(jì)算時(shí)間受DSMC方法中特征參數(shù)(時(shí)間步長(zhǎng)、數(shù)目權(quán)重、網(wǎng)格數(shù)目)的影響規(guī)律,得到以下結(jié)論:
a. 相同聲場(chǎng)作用時(shí)間下,頻率越大,顆粒容積份額變化越迅速,顆粒空間分布也越不均勻;聲場(chǎng)頻率對(duì)碰撞率有顯著影響,而對(duì)計(jì)算時(shí)間的影響很小.
b. 取樣顆粒數(shù)目權(quán)重增加,顆粒碰撞率的變化很小,計(jì)算時(shí)間則迅速降低;網(wǎng)格數(shù)目增加,碰撞率降低,計(jì)算時(shí)間先迅速降低,而后低頻時(shí)基本不變,高頻時(shí)略有上升.
c. 隨著碰撞時(shí)間步長(zhǎng)的增加,低頻聲場(chǎng)中顆粒碰撞率單調(diào)增加,高頻聲場(chǎng)中顆粒碰撞率先增加而后出現(xiàn)波動(dòng);計(jì)算時(shí)間減少,減少量在碰撞步長(zhǎng)較小時(shí)最為明顯.
[1] ZHANG L W,CHEN X,XUE X D,et al.Long-term exposure to high particulate matter pollution and cardiovascular mortality:a 12-year cohort study in four cities in northern China[J].Environment International,2014,62:41-47.
[2] WANG T J,JIANG F,DENG J J,et al.Urban air quality and regional haze weather forecast for Yangtze River Delta region[J].Atmospheric Environment,2012,58:70-83.
[3] WEN C,XU M H,YU D X,et al.PM10formation during the combustion of N2-char and CO2-char of Chinese coals[J].Proceedings of the Combustion Institute,2013,34(2):2383-2392.
[4] EHRLICH C,NOLL G,KALKOFF W D,et al.PM10,PM2.5and PM1.0—Emissions from industrial plants—results from measurement programmes in Germany[J].Atmospheric Environment,2007,41(29):6236-6254.
[5] WANG J,HU Z M,CHEN Y Y,et al.Contamination characteristics and possible sources of PM10and PM2.5in different functional areas of Shanghai,China[J].Atmospheric Environment,2013,68:221-229.
[6] 楊林軍.燃燒源細(xì)顆粒物污染控制技術(shù)[M].北京:化學(xué)工業(yè)出版社,2011.
[7] 胡惠敏,李瑞陽(yáng),蔡萌,等.聲波與其他方法聯(lián)合作用脫除細(xì)顆粒物的研究進(jìn)展[J].上海理工大學(xué)學(xué)報(bào),2016,38(1):13-18.
[8] 張明俊,凡鳳仙.細(xì)顆粒物的聲凝并數(shù)值模擬研究進(jìn)展[J].化工進(jìn)展,2012,31(8):1671-1676.
[9] SHENG C D,SHEN X L.Modelling of acoustic agglomeration processes using the direct simulation Monte Carlo method[J].Journal of Aerosol Science,2006,37(1):16-36.
[10] SHENG C D,SHEN X L.Simulation of acoustic agglomeration processes of poly-disperse solid particles[J].Aerosol Science and Technology,2007,41(1):1-13.
[11] 凡鳳仙,袁竹林.外加聲場(chǎng)對(duì)增加PM2.5碰撞幾率的數(shù)值模擬研究[J].中國(guó)電機(jī)工程學(xué)報(bào),2006,26(11):12-16.
[12] 凡鳳仙,袁竹林,趙兵,等.駐波聲場(chǎng)中細(xì)微顆粒凝并的數(shù)值模擬[J].燃燒科學(xué)與技術(shù),2008,14(3):253-258.
[13] FAN F X,YANG X F,KIM C N.Direct simulation of inhalable particle motion and collision in a standing wave field[J].Journal of Mechanical Science and Technology,2013,27(6):1707-1712.
[14] WEI J M.A fast Monte Carlo method based on an acceptance-rejection scheme for particle coagulation[J].Aerosol and Air Quality Research,2013,13:1273-1281.
[15] HE Y X,ZHAO H B,WANG H M,et al.Differentially weighted direct simulation Monte Carlo method for particle collision in gas-solid flows[J].Particuology,2015,21:135-145.
[16] KARCHANI A,MYONG R S.Convergence analysis of the direct simulation Monte Carlo based on the physical laws of conservation[J].Computers & Fluids,2015,115:98-114.
[17] 宋曉通,凡鳳仙.駐波聲場(chǎng)中可吸入顆粒物漂移的影響因素分析[J].熱能動(dòng)力工程,2016,31(1):81-86.
[18] CLECKLER J,ELGHOBASHI S,LIU F.On the motion of inertial particles by sound waves[J].Physics of Fluids,2012,24(3):033301.
[19] TSUJI Y,TANAKA T,YONEMURA S.Cluster patterns in circulating fluidized beds predicted by numerical simulation(discrete particle model versus two-fluid model)[J].Powder Technology,1998,95(3):254-264.
[20] ROSTAMI M,ARDESHIR A,AHMADI G,et al.On the effect of gravitational and hydrodynamic forces on particle motion in a quiescent fluid at high particle Reynolds numbers[J].Canadian Journal of Physics,2008,86(6):791-799.
(編輯:丁紅藝)
Influence of Characteristic Parameters in the Simulation of Acoustic Particle Collision Using DSMC Method
GUO Xiaowu1,2, FAN Fengxian1,2, HU Xiaohong1,2, SU Mingxu1,2
(1.School of Energy and Power Engineering,University of Shanghai for Science and Technology,Shanghai 200093,China; 2.Shanghai Key Laboratory of Multiphase Flow and Heat Transfer in Power Engineering,University of Shanghai for Science and Technology,Shanghai 200093,China)
An inter-particle collision model was established based on the direct simulation Monte Carlo (DSMC) method as well as the force analysis and the motion equation.The collision process of the particles in acoustic field was examined.Through changing the numerical simulation conditions,the influences of the characteristic parameters,such as the weight value of the number of sampling particles,the cell number and the time-step size used in DSMC method,on the particle collision rate and computational time were discussed.The results show that the higher the acoustic frequency is,the more rapidly the particle volume fraction changes and the more unevenly the particles are distributed in space.The increase of the weight value of the number has a small effect on the collision rate,but it significantly reduces the computational time.As the cell number increases,the collision rate decreases,whereas the computational time decreases rapidly at first and then almost keeps constant at lower acoustic frequency and increases slightly at higher acoustic frequency.It is also found that as the time-step size for determining the inter-particle collision increases,the collision rate increases monotonously at lower frequency,while it increases at first and then fluctuates at higher frequency.The increase of the time-step size results in the decrease of the computational time and the decrease is more obvious with a smaller time-step size.
acoustic field; particle; collision; direct simulation Monte Carlo method; characteristic parameter
1007-6735(2016)05-0419-08
10.13255/j.cnki.jusst.2016.05.003
2016-01-25
國(guó)家自然科學(xué)基金資助項(xiàng)目(51206113,51176128,51576130);上海市科委科研計(jì)劃資助項(xiàng)目(13DZ2260900)
郭孝武(1988-),男,碩士研究生.研究方向:氣固兩相流數(shù)值模擬.E-mail:ross_2015@163.com
凡鳳仙(1982-),女,副教授.研究方向:燃燒源污染物排放控制、氣固流動(dòng)與傳遞.
E-mail:fanfengxian@hotmail.com
X 513
A
上海理工大學(xué)學(xué)報(bào)2016年5期