車 威，魏明銳，李 松

(1.武漢理工大學(xué)現(xiàn)代汽車零部件技術(shù)湖北省重點(diǎn)實(shí)驗(yàn)室，武漢 430070；2.武漢理工大學(xué)汽車零部件技術(shù)湖北省協(xié)同創(chuàng)新中心，武漢 430070； 3.黃淮學(xué)院電子科學(xué)與工程系，駐馬店 463000)

2016003

基于蒙特卡洛方法的碳煙生長(zhǎng)過(guò)程數(shù)值模擬*

車 威1,2,3，魏明銳1,2，李 松1,2

(1.武漢理工大學(xué)現(xiàn)代汽車零部件技術(shù)湖北省重點(diǎn)實(shí)驗(yàn)室，武漢 430070；2.武漢理工大學(xué)汽車零部件技術(shù)湖北省協(xié)同創(chuàng)新中心，武漢 430070； 3.黃淮學(xué)院電子科學(xué)與工程系，駐馬店 463000)

基于彈射凝聚模型構(gòu)建碳煙顆粒間碰撞的物理模型，采用蒙特卡洛方法構(gòu)建碳煙顆粒凝聚和表面生長(zhǎng)的數(shù)學(xué)模型并進(jìn)行數(shù)值模擬。引入分形維數(shù)、包容半徑和凝聚體體積來(lái)描述凝聚體分形生長(zhǎng)特征，分析表面生長(zhǎng)過(guò)程對(duì)這些參量的影響。對(duì)比了不同環(huán)境壓力和溫度對(duì)碰撞頻率、粒徑增長(zhǎng)率和分形維數(shù)的影響。最后通過(guò)枝接過(guò)程的模擬，獲得了鏈狀凝聚體的形貌特征，并與真實(shí)的碳煙照片進(jìn)行對(duì)比，結(jié)果表明，枝接過(guò)程模擬出的碳煙結(jié)構(gòu)更接近于尾氣中真實(shí)的鏈狀碳煙形貌。

碳煙顆粒；彈射凝聚；蒙特卡洛方法；分形維數(shù)

前言

柴油機(jī)工作過(guò)程中產(chǎn)生的廢氣中顆粒物(PM)等對(duì)人類健康和環(huán)境有極大危害。碳煙顆粒的生成要經(jīng)歷燃料熱解、芳香烴形成與生長(zhǎng)、碳核形成和顆粒生長(zhǎng)4個(gè)階段[1]，其間包含有復(fù)雜的物理化學(xué)過(guò)程，涉及從氣相分子到固相顆粒的相變過(guò)程。目前氣相反應(yīng)的計(jì)算可以采用詳細(xì)的化學(xué)反應(yīng)機(jī)理[2]，固體顆粒動(dòng)力學(xué)過(guò)程主要包括成核后顆粒的碰撞、凝聚、團(tuán)聚和破碎等過(guò)程，這些過(guò)程都非常復(fù)雜，加之顆粒動(dòng)力學(xué)過(guò)程在空間上的非連續(xù)性和時(shí)間上的非穩(wěn)定性，難以用確定的數(shù)學(xué)關(guān)系式進(jìn)行表達(dá)[3]，目前主要是采用顆粒群平衡方程來(lái)描述，該方程是一個(gè)積分微分方程，通常沒(méi)有解析解[4]，其求解方法主要有有限元法、有限差分法、矩方法和蒙特卡洛方法等。而針對(duì)碳煙顆粒的成核過(guò)程，目前在實(shí)驗(yàn)上還沒(méi)有很好的手段，只能采用理論與經(jīng)驗(yàn)?zāi)Ｐ湍M顆粒成核?，F(xiàn)有的成核機(jī)理大致可分為3種[5]：類似富勒烯結(jié)構(gòu)的成核機(jī)理、多環(huán)芳香烴(PAHs)通過(guò)物理凝聚的成核機(jī)理和PAHs通過(guò)化學(xué)環(huán)化反應(yīng)的成核機(jī)理。

文獻(xiàn)[6]中雖然用矩方法獲得了碳煙顆粒動(dòng)力學(xué)演變的相關(guān)信息[6]，但是沒(méi)有獲得碳煙顆粒的形貌特征。文獻(xiàn)[7]中采用分形理論對(duì)碳煙的形貌特征進(jìn)行了分析和模擬，但沒(méi)有研究表面生長(zhǎng)和環(huán)境因素對(duì)分形維數(shù)的影響。文獻(xiàn)[8]和文獻(xiàn)[9]中采用蒙特卡洛方法進(jìn)行的計(jì)算集中在求解顆粒群平衡方程，包含顆粒的成核、凝并和氧化過(guò)程，而本文中則采用蒙特卡洛方法來(lái)構(gòu)建碳煙顆粒凝聚和表面生長(zhǎng)的數(shù)學(xué)模型并進(jìn)行數(shù)值求解。

首先利用分形維數(shù)、包容半徑和凝聚體體積等參量分析粒子的凝聚特點(diǎn)，其次探究環(huán)境壓力和環(huán)境溫度對(duì)碰撞頻率、粒徑增長(zhǎng)率和分形維數(shù)的影響，最后模擬兩種情況下的枝接過(guò)程，并與真實(shí)碳煙照片進(jìn)行對(duì)比，探討碳煙顆粒的形成機(jī)制。

1 碳煙顆粒動(dòng)力學(xué)理論與建模

本文中不考慮成核前的氣相反應(yīng)過(guò)程，假設(shè)初始碳核大小相同，并在空間隨機(jī)產(chǎn)生，然后按照彈射運(yùn)動(dòng)形式與其他粒子碰撞，在碰撞過(guò)程中考慮表面反應(yīng)和氧化反應(yīng)。

1.1 初始粒子與凝聚體碰撞——凝聚過(guò)程

粒子的運(yùn)動(dòng)形式有布朗運(yùn)動(dòng)、彈射運(yùn)動(dòng)和受限擴(kuò)散等，本文中選取彈射運(yùn)動(dòng)模式來(lái)模擬粒子間碰撞凝聚。

如圖1(a)所示，將笛卡爾坐標(biāo)系原點(diǎn)定義為凝聚中心，在凝聚中心設(shè)置一個(gè)粒子，定義該粒子為凝聚體，利用蒙特卡洛方法在空間中隨機(jī)產(chǎn)生一個(gè)初始粒子，并以隨機(jī)方向沿直線彈射運(yùn)動(dòng)，與凝聚體之間產(chǎn)生相對(duì)運(yùn)動(dòng)。若發(fā)生碰撞則粘附在一起并作為下次碰撞的凝聚體，如圖1(b)所示；若不發(fā)生碰撞，初始粒子沿著隨機(jī)運(yùn)動(dòng)方向超出定義邊界后消失，然后以相同方式產(chǎn)生下一個(gè)初始粒子，這樣周而復(fù)始，最終完成微粒生長(zhǎng)模擬過(guò)程。初始粒子與凝聚體位置的關(guān)系有一種特殊情況，即初始粒子的行走軌跡恰好與凝聚體相切，在本文中對(duì)這種情況也按照粘附在一起處理。

圖1 初始粒子與凝聚體碰撞過(guò)程

1.2 凝聚體與凝聚體碰撞——枝接過(guò)程

圖2為凝聚體與凝聚體之間碰撞的枝接過(guò)程。當(dāng)發(fā)生枝接碰撞時(shí)，每個(gè)凝聚體所含碳煙粒子數(shù)不再增加，枝接只和兩個(gè)凝聚體尺寸大小和初始凝聚體的運(yùn)動(dòng)軌跡有關(guān)。

圖2 枝接過(guò)程

1.3 粒子半徑增長(zhǎng)量

粒子的表面反應(yīng)包括粒子的表面生長(zhǎng)和表面氧化兩個(gè)過(guò)程。氣相組分在粒子表面的生長(zhǎng)和沉積使得粒子半徑增大，而由OH和O2引起的氧化反應(yīng)則使粒子半徑減小，粒子最終直徑的大小是表面生長(zhǎng)和氧化過(guò)程共同作用的結(jié)果。本文中根據(jù)粒徑增長(zhǎng)率和碰撞時(shí)間計(jì)算出粒子表面生長(zhǎng)量。

Knudsen數(shù)定義為Kn=2λ/d,λ為系統(tǒng)中粒子平均自由程，d為粒子大小的特征尺度，如粒子的直徑。粒子碰撞過(guò)程根據(jù)Kn數(shù)的不同可分為3種模式，分別為

(1)Kn→0：連續(xù)介質(zhì)模式

(2)Kn→：自由分子模式

(3)0.25

連續(xù)介質(zhì)模式下，碰撞頻率[10]為

(1)

其中K=2kBT/(3η)式中：kB為玻爾茲曼常數(shù)；T為熱力學(xué)溫度；η為氣體分子黏度；Ci為連續(xù)介質(zhì)狀態(tài)的修正系數(shù)，Ci=1+1.257Kni；mi為凝聚體i的質(zhì)量；Cj，mj的含義類似。

自由分子模式下，碰撞頻率為

(2)

式中：mc為碳原子質(zhì)量；ρ為碳煙粒子密度。

介于連續(xù)模式與自由分子模式之間的為過(guò)渡模式。通常情況下柴油機(jī)燃燒壓力為6.5MPa，溫度為1 500K，對(duì)于直徑分別為1，2，5和10nm的粒子，Kn數(shù)分別為5.2，2.6，1.0和0.5，其值介于0.25和10之間[11]，為過(guò)渡模式。過(guò)渡模式碰撞頻率可取自由模式和連續(xù)模式碰撞頻率的平均值，無(wú)量綱化后的表達(dá)式[7]為

(3)

其中

第i個(gè)和第j個(gè)凝聚體發(fā)生第n次碰撞所經(jīng)歷的時(shí)間為

(4)

式中：N為第n次碰撞時(shí)刻的粒子數(shù)密度。凝聚經(jīng)歷的總時(shí)間為tn+1=tn+Δti,j。

粒子的表面生長(zhǎng)反應(yīng)促進(jìn)了粒子半徑的增加，進(jìn)而使整個(gè)凝聚體半徑增加，粒子半徑R隨時(shí)間t的變化率可用粒徑增長(zhǎng)率來(lái)描述，即

dR/dt=Ω

(5)

式中：Ω為粒徑增長(zhǎng)率。黑體字符代表有量綱的變量，下同。將式(5)的各個(gè)參量無(wú)量綱化處理，令

(6)

式中：R0為出現(xiàn)在模型中第一個(gè)粒子的半徑；t0為粒子間發(fā)生第一次碰撞經(jīng)歷的時(shí)間。對(duì)式(5)進(jìn)行變換：

(7)

結(jié)合式(5)和式(6)，式(7)可變?yōu)?/p>

(8)

結(jié)合式(8)，R對(duì)于t的1階泰勒展開(kāi)式為

即R(t+Δt)-R(t)=ΩΔt，或?qū)憺?/p>

ΔR=ΩΔt

(9)

ΔR即為每次碰撞事件發(fā)生時(shí)粒子表面增加的厚度。

2 凝聚體參量的蒙特卡洛算法

由于粒子表面增長(zhǎng)，碳煙粒子間存在交叉重疊，在計(jì)算凝聚體的體積、表面積和回轉(zhuǎn)半徑等參量時(shí)，如果只是把每個(gè)粒子對(duì)應(yīng)參量值簡(jiǎn)單相加，則隨著碳煙粒子數(shù)的增加，計(jì)算誤差逐漸增大。為了減小這種計(jì)算誤差，采用蒙特卡洛取樣算法[12]。定義碳煙粒子數(shù)為n的凝聚體中每個(gè)碳煙粒子為B1，B2，…，Bn，隨機(jī)選擇凝聚體內(nèi)的一個(gè)碳煙粒子Bn，并在此粒子內(nèi)部隨機(jī)選擇一個(gè)點(diǎn)pm，建立一個(gè)關(guān)于pm,Bn的函數(shù)f(pm,Bn)，m為粒子Bn內(nèi)隨機(jī)選擇的點(diǎn)的序號(hào)，根據(jù)粒子數(shù)n等信息計(jì)算出取樣點(diǎn)的數(shù)量M。當(dāng)pm點(diǎn)落入Bn編號(hào)為最小的粒子內(nèi)時(shí)，此時(shí)f(pm,Bn)=1，當(dāng)pm點(diǎn)落入Bn編號(hào)不是最小的粒子內(nèi)時(shí)，則f(pm,Bn)=0，粒子就不被計(jì)數(shù)，即

f(pm,Bn)=

(10)

凝聚體表面積計(jì)算采用類似的算法，不同之處在于所選擇的點(diǎn)必須位于粒子表面上。

(11)

則凝聚體回轉(zhuǎn)半徑Rg為

(12)

3 分形維數(shù)

歐氏空間用1，2，3這樣的整數(shù)表示線、面、空間的維數(shù)，卻無(wú)法表示非規(guī)則、自相似的幾何圖形的維數(shù)，比如測(cè)定海岸線的長(zhǎng)度等。分形維數(shù)可以表示分形集的不規(guī)則程度，突破了維數(shù)為整數(shù)的界限，一般為分?jǐn)?shù)。

隨著研究的深入，分形理論逐漸應(yīng)用到物理學(xué)、化學(xué)、生物學(xué)、天文學(xué)和經(jīng)濟(jì)學(xué)等，為許多問(wèn)題的解決提供了新的思路和方法。分形維數(shù)的表示方法有很多，本文中采用的分形維數(shù)Df表示為

(13)

Df=slope{lnri, lnVri}

(14)

4 計(jì)算結(jié)果與分析

本文中應(yīng)用C++語(yǔ)言編程構(gòu)建了粒子碰撞凝聚的計(jì)算模型。以柴油機(jī)工作過(guò)程的溫度、壓力、粒子數(shù)密度、C2H2沉積速率、O2和OH的氧化速率等作為初始數(shù)據(jù)[7]，對(duì)碳煙粒子碰撞凝聚生長(zhǎng)過(guò)程進(jìn)行了數(shù)值模擬。

4.1 凝聚碰撞

在環(huán)境壓力保持不變的情況下，分3種情況對(duì)比表面反應(yīng)對(duì)粒子凝聚過(guò)程的影響：C1無(wú)表面生長(zhǎng)和氧化；C2只有表面生長(zhǎng)；C3有表面生長(zhǎng)和氧化。

圖3是在粒子數(shù)為400時(shí)計(jì)算得到的上述3種情況下凝聚體形態(tài)圖。通過(guò)對(duì)比可知，對(duì)于圖3(a)，粒子直徑始終不發(fā)生變化，凝聚體內(nèi)部粒子大小和外部粒子大小一致；對(duì)于圖3(b)，只存在表面生長(zhǎng)，凝聚體內(nèi)部粒子先于外部粒子粘附在凝聚體上，隨著碰撞過(guò)程進(jìn)行，內(nèi)部粒子在空間環(huán)境中存在的時(shí)間較長(zhǎng)，粒徑增長(zhǎng)較大，所以粒子交叉重疊現(xiàn)象明顯；而外部粒子在空間環(huán)境中存在時(shí)間較短，所以粒徑增加較小，交叉重疊現(xiàn)象較弱。對(duì)于圖3(c)，同時(shí)存在表面生長(zhǎng)和氧化反應(yīng)，凝聚體交叉重疊情況相對(duì)C2較弱，凝聚體包容半徑相對(duì)于C2也較小，但由于氧化作用在整個(gè)粒子凝聚過(guò)程中所占比重較小，與圖3(b)相比不是特別明顯，這在圖4中進(jìn)一步得到了證實(shí)。

圖3 3種情況下凝聚體的形態(tài)

圖4是在粒子數(shù)n為750，環(huán)境壓力為1MPa的情況下凝聚體參量的變化，由于存在隨機(jī)波動(dòng)，故進(jìn)行了50次重復(fù)計(jì)算，然后對(duì)計(jì)算結(jié)果進(jìn)行平均處理。通過(guò)圖4(a)看出，C1與C2和C3情況相比，凝聚體分形維數(shù)較小。這是因?yàn)?，一方面由彈射凝聚模型特點(diǎn)造成凝聚體結(jié)構(gòu)較松散，另一方面粒子間不存在交叉重疊情況；另外C1曲線后期略有上升是因?yàn)橛猩倭苛Ｗ舆M(jìn)入凝聚體內(nèi)部導(dǎo)致凝聚體變得更加致密。在C2和C3兩種情況下，在碰撞初期分形維數(shù)變化較大，主要是因?yàn)楫?dāng)?shù)谝粋€(gè)初始粒子出現(xiàn)在空間成為凝聚體時(shí)，分形維數(shù)為3。但當(dāng)少數(shù)粒子碰撞粘附在一起后，凝聚體結(jié)構(gòu)較松散，粒子間交叉重疊不明顯，分形維數(shù)較低。而隨著粒子表面生長(zhǎng)和氧化反應(yīng)的繼續(xù)進(jìn)行和粒子的不斷填充，凝聚體更加致密和接近于團(tuán)簇狀，使分形維數(shù)又有較大提高。C3與C2相比，還存在粒子的表面氧化，使粒子的密實(shí)程度降低，分形維數(shù)也有一定程度降低。圖4(b)是凝聚體的最小包容半徑隨粒子數(shù)的變化，可知，最小包容半徑初期增加較快，后期增加較慢，主要是初期凝聚體半徑相對(duì)初始粒子半徑比例較小，所以粘附到凝聚體的初始粒子體積對(duì)凝聚體體積影響較明顯，當(dāng)凝聚體粒子數(shù)較多時(shí)，整個(gè)凝聚體包容半徑達(dá)到30以上，而此時(shí)初始粒子半徑只有1，對(duì)凝聚體包容半徑影響較小，所以隨著粒子繼續(xù)增加，凝聚體包容半徑增長(zhǎng)減慢。圖4(c)為凝聚體體積隨粒子數(shù)的變化，對(duì)于C1情況，凝聚體體積就是單個(gè)粒子體積乘以粒子數(shù)，所以圖線呈線性增長(zhǎng)；而對(duì)于C2和C3，凝聚體體積的增大除了粒子增長(zhǎng)的因素之外，表面生長(zhǎng)的影響也較大，這可以從圖4中看出，所以C2和C3曲線與C1曲線的差別較大；而C2和C3曲線之間也有差別，說(shuō)明表面氧化對(duì)凝聚體體積也有影響。由于C2和C3兩種情況后期表面反應(yīng)都大大減弱，粒子半徑增大越來(lái)越不明顯，所以兩條曲線后期斜率基本保持不變，且凝聚體體積增長(zhǎng)和粒子數(shù)近似呈線性關(guān)系。

圖4 凝聚體參量值在3種情況下的對(duì)比

圖5 環(huán)境壓力變化對(duì)凝聚過(guò)程參量的影響

圖6為粒子碰撞過(guò)程中有表面生長(zhǎng)和氧化，環(huán)境壓力保持在1MPa，不同溫度下粒子的凝聚過(guò)程。由圖6(a)可見(jiàn)，隨著溫度增加，粒子碰撞頻率減小。由圖6(b)可見(jiàn)，由于C2H2沉積速率先增大后減小最后趨于穩(wěn)定，O2和OH的氧化速率一直減小最后趨于穩(wěn)定，且表面生長(zhǎng)速率大于表面氧化速率，兩者共同作用的結(jié)果使得粒徑增長(zhǎng)率先增大后減小，對(duì)應(yīng)圖上出現(xiàn)了一個(gè)峰值，后期表面生長(zhǎng)速率趨于穩(wěn)定而氧化速率的影響較弱，且后期碰撞粒子直徑較小，所以對(duì)應(yīng)粒徑增長(zhǎng)率緩慢增加。溫度升高導(dǎo)致碳煙表面氧化速率比生長(zhǎng)速率增加的更快，故粒徑增長(zhǎng)率減小。圖6(c)為凝聚體分形維數(shù)隨粒子數(shù)變化關(guān)系，整體變化趨勢(shì)同圖5(c)。由于粒徑增長(zhǎng)率隨溫度升高而降低，致使凝聚體致密程度降低，所以分形維數(shù)減小。

圖6 溫度變化對(duì)凝聚過(guò)程參量的影響

4.2 枝接碰撞

圖7是通過(guò)枝接形成的凝聚體的形態(tài)圖，發(fā)生枝接凝聚的4個(gè)小凝聚體所含粒子數(shù)分別為106，119，121和137，總和為483個(gè)。圖7(a)為無(wú)表面生長(zhǎng)和氧化過(guò)程，圖7(b)為伴隨有表面生長(zhǎng)和氧化過(guò)程。可以看出，通過(guò)枝接過(guò)程形成的凝聚體具有明顯的鏈狀結(jié)構(gòu)，同樣，表面生長(zhǎng)和氧化可以改變凝聚體的致密程度和分形維數(shù)等參數(shù)，表1給出了凝聚體的相應(yīng)參數(shù)。

圖7 兩種情況下枝接形成的凝聚體形貌

分類粒子數(shù)n分形維數(shù)Df包容半徑Re回轉(zhuǎn)半徑Rg無(wú)表面生長(zhǎng)和氧化4832.02432.87216.790有表面生長(zhǎng)和氧化4832.25066.80327.031

圖8為柴油機(jī)排放尾氣中的碳煙顆粒物形貌照片[13]，其中大部分為鏈狀結(jié)構(gòu)，而且較大的鏈狀結(jié)構(gòu)中包含著較小的鏈狀結(jié)構(gòu)，較小的鏈狀結(jié)構(gòu)中還包含有更小的鏈狀結(jié)構(gòu)，呈現(xiàn)分形的自相似性。由圖3、圖7和圖8對(duì)比可知，通過(guò)枝接過(guò)程模擬出的碳煙結(jié)構(gòu)更接近于尾氣中真實(shí)的鏈狀碳煙形貌并具有分形特征。

圖8 真實(shí)柴油機(jī)排放產(chǎn)生的碳煙形貌

5 結(jié)論

應(yīng)用彈射凝聚模型和蒙特卡洛方法構(gòu)建了碳煙顆粒凝聚和表面生長(zhǎng)的數(shù)學(xué)物理模型，并以分形維數(shù)描述凝聚體生長(zhǎng)過(guò)程和形貌特征，根據(jù)不同條件進(jìn)行一系列數(shù)值求解，得到以下結(jié)論。

(1) 碳煙的表面生長(zhǎng)會(huì)明顯提高凝聚體的分形維數(shù)、包容半徑和凝聚體體積，并使凝聚體變得更加致密，表面氧化則起相反作用。

(2) 隨著環(huán)境壓力的升高，粒子碰撞頻率減小且波動(dòng)性減弱，導(dǎo)致分形維數(shù)增大；隨著環(huán)境溫度的升高，粒子的粒徑增長(zhǎng)率減小，導(dǎo)致分形維數(shù)降低。

(3) 通過(guò)與柴油機(jī)尾氣中真實(shí)碳煙照片對(duì)比可知，枝接過(guò)程模擬出的碳煙結(jié)構(gòu)更接近于尾氣中真實(shí)的鏈狀碳煙形貌并具有分形特征。

通過(guò)上述結(jié)論可以看出，利用碳煙顆粒的形貌特征能夠獲得燃燒環(huán)境和燃燒過(guò)程的相關(guān)信息，如表面反應(yīng)、壓力和溫度影響等情況，從而為優(yōu)化燃燒、降低碳煙提供依據(jù)。

[1] FENG T, REITZ R D, FOSTER D E, et al. Nine-step Phenomenological Diesel Soot Model Validated over a Wide Range of Engine Conditions[J]. International Journal of Thermal Sciences, 2009, 48(6):1223-1234.

[2] 陳亮, 成曉北, 顏方沁,等. 基于改進(jìn)的詳細(xì)碳煙模型的柴油燃燒碳煙顆粒物的生成特性[J].燃燒科學(xué)與技術(shù), 2013, 19(3): 234-240.

[3] 趙海波, 鄭楚光. 離散系統(tǒng)動(dòng)力學(xué)演變過(guò)程的顆粒群平衡模擬[M]. 北京: 科學(xué)出版社, 2008.

[4] PATTERSON R I A, WAGNER W, KRAFT M. Stochastic Weighted Particle Methods for Population Balance Equations[J]. Journal of Computational Physics, 2011, 230(19): 7456-7472.

[5] WANG H. Formation of Nascent Soot and Other Condensed-phase Materials in Flames[J]. Proceedings of the Combustion Institute, 2011, 33(1): 41-67.

[6] FRENKLACH M, HARRIS S J. Aerosol Dynamics Modeling Using the Method of Moments[J]. Journal of Colloid and Interface Science, 1987, 118(1):252-261.

[7] MITCHELL P A. Monte Carlo Simulation of Soot Aggregation with Simultaneous Surface Growth[D]. University of California, Berkeley, 2001.

[8] CELNIK M, RAJ A, WEST R, et al.Aromatic Site Description of Soot Particles[J]. Combustion and Flame, 2008,155(1-2):161-180.

[9] RAJ A, SANDER M, JANARDHANAN V, et al.A Study on the Coagulation of Polycyclic Aromatic Hydrocarbon Clusters to Determine Their Collision Efficiency[J]. Combustion and Flame, 2010,157(3):523-534.

[10] KAZAKOV A, FRENKLACH M. Dynamic Modeling of Soot Particle Coagulation and Aggregation: Implementation with the Method of Moments and Application to High-pressure Laminar Premixed Flames[J]. Combustion and Flame, 1998, 114(3-4): 484-501.

[11] 崔心存, 金國(guó)棟. 內(nèi)燃機(jī)排放凈化[M]. 武漢：華中理工大學(xué)出版社，1990.

[12] KARP R M, LUBY M, MADRAS N.Monte Carlo Approximation Algorithms for Enumeration Problems[J]. Journal of Algorithms in Cognition, Informatics and Logic, 1989, 10(3):429-448.

[13] Van GULIJK C, MARIJNISSEN J C M, MAKKEE M, et al. Measuring Diesel Soot with a Scanning Mobility Particle Sizer and an Electrical Low-pressure Impactor: Performance Assessment with a Model for Fractal-like Agglomerates[J]. Journal of Aerosol Science, 2004, 35(5):633-655.

Numerical Simulation of Soot Growth Process Based on Monte Carlo Method

Che Wei1,2,3, Wei Mingrui1,2& Li Song1,2

1.HubeiKeyLaboratoryofAdvancedTechnologyforAutomotiveComponents,WuhanUniversityofTechnology,Wuhan430070;2.HubeiCollaborativeInnovationCenterforAutomotiveComponentsTechnology,WuhanUniversityofTechnology,Wuhan430070;3.DepartmentofElectronicsandScienceEngineering,HuanghuaiUniversity,Zhumadian463000

The physical model for the collision between soot particles is established based on ballistic aggregation model, and Monte Carlo method is used to build the mathematical model for soot particle aggregation and surface growth, with numerical simulation conducted. Fractal dimension, enclosing radius and aggregate volume are introduced to describe the fractal growth characteristics of aggregates, and the effects of surface growth process on these parameters are analyzed. The effects of different ambient pressures and temperatures on the collision frequency, the growth rate of particle size, and fractal dimension are compared. Finally through a simulation on coalescence process, the morphological feature of chain-like aggregates is obtained and compared with real soot photos. The results show that the soot structure simulated by coalescence process is closer to real morphology of chain-like soot particles in exhaust gas.

soot particles; ballistic aggregation; Monte Carlo method; fractal dimension

*國(guó)家自然科學(xué)基金(51276132)和中央高校基本科研業(yè)務(wù)費(fèi)專項(xiàng)資金(2015Ⅲ040)資助。

原稿收到日期為2015年5月14日，修改稿收到日期為2015年8月4日。