煙塵氣溶膠光散射特性的數(shù)值模擬研究

侯宏錄，劉 凱,趙群英

(西安工業(yè)大學 光電工程學院，西安 710021)

?

煙塵氣溶膠光散射特性的數(shù)值模擬研究

侯宏錄，劉 凱,趙群英

(西安工業(yè)大學 光電工程學院，西安 710021)

為研究光在煙塵氣溶膠介質(zhì)中的傳輸特性，通過Mie散射算法對單個煙塵顆粒的光學參數(shù)(單次反照率、不對稱因子等)進行了數(shù)值計算.基于輻射輸運理論，采用蒙特卡羅法模擬了光在煙塵氣溶膠介質(zhì)中的傳輸過程，分析研究了在不同顆粒粒徑、入射波長以及粒子數(shù)密度條件下，入射光通過煙塵氣溶膠介質(zhì)的透射率隨光學厚度以及散射光強隨散射角的變化規(guī)律.研究結(jié)果表明：入射波長增大，煙塵顆粒的散射能力減弱，通過煙塵氣溶膠介質(zhì)的光強透射率增大；煙塵顆粒粒徑增大，煙塵顆粒的散射能力增強，通過煙塵氣溶膠介質(zhì)的光強透射率減小，且當煙塵顆粒粒徑較大(如1 μm)時，煙塵顆粒的光學特性變化平緩.

煙塵氣溶膠；Mie散射；蒙特卡羅；前向散射；

近年來環(huán)境污染嚴重影響著人類的健康，煙塵氣溶膠物質(zhì)是環(huán)境污染的主要來源.氣溶膠是指漂浮在空氣中，由液體和固體微粒與氣體載體共同組成的多相體系[1].對煙塵氣溶膠排放的檢測是掌握煙塵污染情況以及預測煙塵發(fā)展趨勢的重要手段.由于煙塵氣溶膠的分布存在很大的不確定性[2]，到目前為止，人們對氣溶膠本身的形成和發(fā)展過程的了解還不夠深入，很難直接對氣溶膠各種特性進行直接觀測.所以對煙塵氣溶膠排放的檢測，主要是通過研究煙塵微粒物的物理性質(zhì)(例如電學、光學等)得到.

近年來，國內(nèi)外學者在微粒的光學特性研究方面取得了許多成果.文獻[3]成為粒子散射研究領域的經(jīng)典教科書.文獻[4]通過模擬煙塵等特殊環(huán)境，研究了入射波長為0.53 μm、1.06 μm的激光在煙霧環(huán)境中的傳輸特性，得到波長為1.06 μm的激光光源在煙霧中具有更好的傳輸特性的結(jié)果.文獻[5]利用離散坐標法和Monte Carlo法對兩平行層不同粒子系的散射場進行了計算，結(jié)果表明該方法計算結(jié)果更加精確，計算時間也被大大縮減.文獻[6]分析了上海地區(qū)霾期間低層0～2 km高度氣溶膠微物理屬性與地面相對濕度的關系.文獻[7]分析了東中國海域上空氣溶膠的光學特性，結(jié)果表明在該海域氣溶膠大部分為弱吸收類型(約占60%)，同時也存在強散射性類型(約占20%)和強吸收性類型(約占20%).文獻[8]通過對測得的消光系數(shù)[9]和氣溶膠的質(zhì)量濃度建立關聯(lián)，進而通過消光系數(shù)的垂直分布反演氣溶膠質(zhì)量濃度的垂直分布，建立了兩者的相關性聯(lián)系.文獻[10]模擬研究了氣溶膠的光學特性，如單次散射反照率、不對稱因子等光學參數(shù).在粒子散射理論的研究方面，由于我國技術水平和研究經(jīng)費的限制，除了在部分實驗研究中取得一定的研究成果外，研究還不夠深入，主要體現(xiàn)在缺乏對粒子散射場特性和對復雜粒子體系的研究等.

本文基于Mie散射理論及其算法，借助概率抽樣的方法，對單一煙塵顆粒系的散射光強角分布進行計算，分析了不同煙塵顆粒粒徑、入射波長以及粒子數(shù)密度對散射光強的影響，為進一步采用蒙特卡羅方法數(shù)值模擬服從多分散系分布的煙塵氣溶膠粒子的光學特性研究奠定基礎.

1 概率抽樣模型

由激光光源發(fā)出單色光波，將光波處理成光子，對于單一煙塵氣溶膠顆粒系，采用蒙特卡羅方法解決光在煙塵顆粒系中的輻射輸運問題，假定光子與隨機分布粒子的相互作用為彈性散射.光子在隨機介質(zhì)中或被粒子散射或粒子吸收或從介質(zhì)中逃逸出去.其各部分的數(shù)量關系與顆粒系的總吸收截面和總散射截面直接相關，被散射的光在以煙塵顆粒為中心的各個方向上呈不均勻分布，散射光強分布的形狀與煙塵顆粒的大小、折射率有關，可通過Mie散射算法計算得到.光子在隨機游走過程中兩次與煙塵顆粒相碰撞間的行程是顆粒系濃度、顆粒的總消光截面的函數(shù)，每個光子遭受多重散射時，每一次的散射只與前一次散射有關，即光子歷史的狀態(tài)序列可描述為馬爾科夫過程[11].

假定在密閉容器中，煙塵顆粒系在空間上的分布是均勻的，不考慮容器腔體對光子的吸收、反射及散射.在該假定條件下得到光子隨機游走過程的概率抽樣模型，首先設定光子的初始參數(shù)，對光子與煙塵粒子的每次碰撞后的散射路徑和散射方向進行計算，追蹤光子的移動軌跡，直到光子逃逸出煙塵顆粒系或被煙塵粒子吸收后終止追蹤，最終統(tǒng)計每個光子的信息，計算得到散射光強.

1) 初始值設定

初始值設定是指確定光子的初始狀態(tài)，包括光子的位置狀態(tài)、入射角度以及權重因子等.設定光子位置在坐標系原點，初始方向為μ0=cosθ0=1，光子數(shù)目為10 000.

2) 粒子散射路徑抽樣

光子與粒子發(fā)生碰撞而散射，經(jīng)過距離為L的位置上存在的概率為exp(-ctL)，粒子間碰撞距離L的抽樣函數(shù)為

L=-lnξ1/ct

(1)

式中：ξ1為服從(0,1)均勻分布的隨機數(shù);ct為單個粒子的衰減截面和粒子數(shù)密度的乘積.

3)光子散射方向和位置的跟蹤

散射相位函數(shù)取Henyey-Greenstein[12]公式，其函數(shù)形式為

(2)

式中：g為煙塵顆粒的不對稱因子；μ為光子散射方向與入射方向的夾角.

光子的散射角v的抽樣表達式為

(3)

其中ξ2為服從(0,1)均勻分布的隨機數(shù).

考慮各向同性散射，其散射方位角在0～2π均勻分布，因此方位角ψ的抽樣表達式為

ψ=2πξ3

(4)

其中ξ3為服從(0,1)均勻分布的隨機數(shù).

圖1 光子散射方向和微粒位置

圖2 微粒坐標系和原坐標系的散射方向關系

設在xyz系中rm和rm+1的坐標系分別為(x′,y′,z′)和(x,y,z)，建立兩次散射的坐標關系，有

(5)

在原坐標系中方向余弦(U,V,W)為

(6)

(7)

結(jié)合式(1)～(7)，可計算跟蹤光子的散射方向和位置.

4) 光子從介質(zhì)中逃逸或被吸收的檢驗

當光子沿入射光光軸方向移動的距離z大于光學厚度τ或距離z為負值時，表示光子從煙塵介質(zhì)中逃逸出去.如果光子不被煙塵微粒吸收的概率小于光子臨界生存值ωc(取ωc=10-10)，則表示光子被微粒吸收，終止對光子的追蹤.如果光子既沒有被微粒吸收，也沒有逃逸出，則需重新計算煙塵微粒的散射路徑抽樣以及光子位置和散射方向的跟蹤，循環(huán)計算直到光子被煙塵微粒吸收或逃逸出后停止計算.最終統(tǒng)計出射光子的信息,計算得到光強透射率隨煙塵氣溶膠顆粒群厚度及散射光強隨散射角的變化關系.

2 單個顆粒的散射特性

在用抽樣的方法確定光子與煙塵顆粒碰撞后的散射方向時，需要知道單個顆粒的散射相函數(shù)和不對稱因子；確定顆粒的自由行程以及是否被吸收時，需要計算單個顆粒的吸收截面、散射截面、消光截面和單次散射反照率，這些單個顆粒的散射特性可通過Mie散射算法進行求解.

根據(jù)Mie理論，對于單一分散系，即均勻大小的顆粒系，球形粒子的散射相函數(shù)[13]可表示為

(8)

式中：θ為散射角；D為顆粒粒徑；λ為入射波長；α為煙塵顆粒尺寸參數(shù)(α=πD/λ)；Qsca為散射系數(shù)；S1和S2為由貝塞爾函數(shù)和勒讓德函數(shù)組成的無窮級數(shù)，其表達式為

(9)

(10)

式中：x和m分別為煙塵粒子的尺寸參數(shù)和復折射率；μ和μ1分別為周圍介質(zhì)及煙塵顆粒的磁導率.引入兩個函數(shù)：第一類Bessel函數(shù)Jn和第二類Bessel函數(shù)Yn，則jn(z)和hn(1)(z)分別表示為

散射系數(shù)(Qsca)和吸收系數(shù)(Qabs)、消光系數(shù)(Qext)分別為

(11)

單次散射反照率為

ω=Qsca/Qext=1-Qabs/Qext

(12)

單次散射反照率是用來衡量粒子散射和吸收的相對重要性的量，ω越大，表示煙塵粒子的散射能力越強.不對稱因子為

(13)

不對稱因子反映粒子前、后向散射的相對強度，在-1(完全后向散射)到+1(完全前向散射)之間變化.

3 計算結(jié)果和分析

根據(jù)表1不同入射波長情況下煙塵氣溶膠折射指數(shù)的數(shù)據(jù)[14]，計算不同粒徑的煙塵顆粒的不對稱因子、散射系數(shù)、吸收系數(shù)、消光系數(shù)和單次散射反照率，并進行概率抽樣及統(tǒng)計實驗，得到在不同波長、煙塵顆粒粒徑和顆粒粒子數(shù)密度條件下，歸一化散射光強的分布，抽樣實驗的主要結(jié)果.

圖3數(shù)值模擬了在不同顆粒粒徑條件下，單次散射反照率隨波長的變化，如圖3(a)所示，以及在不同波長條件下單次散射反照率隨粒徑的變化，如圖3(b)所示.由圖3可得單次散射反照率隨著入射波長的增加逐漸減??；顆粒粒徑越大，煙塵氣溶膠的單次散射反照率越大；當粒徑大于1 μm左右時，單次反照率的變化趨于平緩.隨著波長的增大，煙塵顆粒的散射能力減弱，吸收能力增強；隨著顆粒粒徑的增大，煙塵顆粒的散射能力增強，吸收能力減弱；當顆粒粒徑較大(比如1 μm)時，煙塵顆粒的散射能力和吸收能力變化平緩.

表1 煙塵氣溶膠折射指數(shù)

圖4模擬了在不同顆粒粒徑條件下，不對稱因子隨波長的變化以及不同波長情況下不對稱因子隨粒徑的變化.由圖4可知,不對稱因子隨入射波長的增加逐漸減?。活w粒粒徑越大，不對稱因子也就越大；當粒徑大于1 μm左右時，單次反照率的變化逐漸不明顯.分析得到,隨著波長的增大，前向散射能力減弱；隨著煙塵顆粒粒徑增大，前向散射能力增強；當煙塵顆粒粒徑較大(比如1 μm)時，前向散射能力基本保持不變.

通過數(shù)字仿真，隨著光學厚度的增大，入射光強的透過率逐漸減少，如圖5所示.由圖5(a)可得,在顆粒粒徑和粒子數(shù)密度不變的情況下，隨著波長λ增大，入射光強的透射率的衰減逐漸減弱；由圖5(b)可得,在波長和粒子數(shù)密度不變的情況下，隨著顆粒粒徑的增大，入射光強的透射率的衰減逐漸增強；由圖5(c)可得,在波長和顆粒粒徑不變的情況下，粒子數(shù)密度越大，入射光強的透射率的衰減越強.

在不同波長、顆粒粒徑和粒子數(shù)密度的情況下，光學厚度為1 m處的歸一化后的散射光強隨散射角的變化關系，如圖6所示，由圖6可知，散射光強的能量主要集中在前向小角度方向上，經(jīng)過統(tǒng)計計算得到前向散射角30°內(nèi)集中了散射光強大約65%的能量.

圖3 單次反照率隨物理參數(shù)的變化

圖4 不對稱因子隨物理參數(shù)的變化

圖5 透射率隨光學厚度的變化

圖6 散射光強隨散射角的變化圖

4 結(jié) 論

通過蒙特卡羅方法在Matlab開發(fā)環(huán)境中數(shù)字模擬了煙塵粒子對光的輸運過程，研究了對于單一煙塵顆粒系.單個顆粒的散射特性和入射光強透射率關于波長、顆粒粒徑以及粒子數(shù)密度三個因素的關系,得出以下結(jié)論：

1) 煙塵的散射能力和前向散射隨著波長的增加減弱，隨著顆粒粒徑的增大而增強；當顆粒粒徑增大到一定程度(如1 μm)時，煙塵顆粒的散射能力變化不明顯.

2) 在同一光學厚度處，入射光強透射率隨入射波長的增大而增大，隨煙塵顆粒粒徑和粒子數(shù)密度的增大而減小.

3) 散射光強的能量主要集中在前向散射的小角度方向上，且前向散射角30°內(nèi)集中了散射光強大約65%的能量.

[1] 貴志成,鄭益群,錢貞成,等.中國東部氣溶膠活化顆粒物濃度對區(qū)域氣候影響的模擬研究[J].南京大學學報(自然科學),2014,50(6):781.

GUI Zhicheng,ZHEN Yiqun,QIAN Zhencheng,et al.Simulation of the Effect on Regional Climate Caused by Aerosol Particle Concentration over Eastern China[J].Journal of Nanjing University(Natural Science),2014,50(6):781.(in Chinese)

[2] 白建輝,王庚辰.黑碳氣溶膠研究新進展[J].科學技術與工程，2005,5(9):585.

BAI Jianhui,WANG Gengchen.New Progress in the Study of Black Carbon Aerocol[J].Science Technology and Engineering,2005,5(9):585.(in Chinese)

[3] MISHCHENKO M I,TRAVIS L D, LACIS A A.Scattering,Absorption,and Emission of Light by Small Particles[M].Cambridge University Press,2002.

[4] 許祖兵.激光大氣傳輸特性分析研究[D].南京:南京理工大學,2006.

XU Zubing.Analysis of Atmospheric Transmission Properties with Laser[D].Nanjing:Nangjing Science and Technology University,2006. (in Chinese)

[5] SOMMERSTEN E R,LOTSBERG J K,STAMNES K,et al.Discrete Ordinate and Monte Carlo Simulations for Polarized Radiative Transfer in a Coupled System Consisting of Two Media with Different Refractive Indices[J].Journal of Quantitative Spectroscopy & Radiative Transfer,2010,111(4):616.

[6] 黃家敏,陳星宇,安靜宇,等.霾期間上海低層氣溶膠微物理屬性與地面相對濕度分析[J].沙漠與綠洲氣象,2015,9(1):39.

HYANG Jiaming,CHEN Xingyu,AN Jingyu,et al.The Analysis of the Low Layer Aerosol Micro Physical Properties and Relative Humidity of Ground during Haze in Shanghai [J].Desert and Oasis Meteorology,2015,9(1):39.(in Chinese)

[7] 王浩然,陳文忠.中國東部海域大氣氣溶膠光學性質(zhì)及其氣溶膠模型分析[J].大氣與環(huán)境光學學報,2013,8(3):179.

WANG Haoran,CHENG Wenzhong.The Analysis of Optical Properties and Model of Aerosol over the East China Sea[J].Journal of Atmospheric and Environmental Optics,2013,8(3):179.(in Chinese)

[8] 韓道文,劉文清,劉建國,等.氣溶膠質(zhì)量濃度空間垂直分布的反演方法[J].2006,33(11):1567.

HAN Daowen,LIU Wenqing,LIU Jianguo,et al.Retrieval Method for Aerosol Mass Concentration Vertical Distribution[J].2006,33(11):1567.(in Chinese)

[9] 王亞偉.光散射理論及其應用技術[M].北京:科學出版社,2013.

WANG Yawei.Light Scattering Theory and Its Application Technology[M].Beijing:Science Press,2013.(in Chinese)

[10] MAJA K,BOX M A,BEAT S,et al.Case Study of Modeled Aerosol Optical Properties During the SAFARI 2000 Campaign[J].Applied Optics,2007,46(22):5263.

[11] 倪麗明.蒙特卡羅方法在煤煙介質(zhì)光分布特性中的應用[D].西安：西安電子科技大學,2009.

NI Liming.Applications of the Monte Carlo Method in the Study of Distributing Characteristics of Light in Soot[D].Xi’an：Xidian University，2009.(in Chinese)

[12] BINZONI T,LEUNG T S,GRANDJBAKHCHE A H,et al.The Use of the Henyey-Greenstein Phase Function in Monte Carlo Simulation in Biomedical Optics[J].Phys Med Biol,2006,51(17):313.

[13] 袁易君,任德明,胡孝勇.Mie理論遞推公式計算散射相位函數(shù)[J].光散射學報，2005，17(4):366.

YUAN Yijun,REN Deming,HU Xiaoyong.Computing Scattering Phase Function by Recursive Formula of Mie Theory[J].Journal of Light Scattering,2005,17(4):366.(in Chinese)

[14] LEVONI C,CERVINO M,GUZZI R,et al.Atmospheric Aerosol Optical Properties:A Database of Radiative Characteristics for Different Components and Classes[J].Applied Optics,1997,36(30):8031

(責任編輯、校對 潘秋岑)

Numerical Simulation of Light-Scattering Property of Soot Aerosol

HOUHonglu,LIUKai,ZHAOQunying

(School of Optoelectronic Engineering,Xi’an Technological University,Xi’an 710021,China)

In order to study the transmission characteristics of light in soot aerosol medium,by scattering algorithm,the single soot particle optical parameters(single albedo,asymmetry factor,etc.) are calculated numerically.Based on the radiation tansport theory,the monte carlo method is used to simulate the light in the process of transmission in soot aerosol medium.The changing of incident light through the transmittance of soot aerosol medium and the varying of scattering intensity with the scattering angle are analyzed,under conditions of different particle size ,incident wavelength and the particle number density.The results show:With the increase of the incidence wavelength,soot particle scattering ability is abate,and the intensity transmissivity through the soot aerosol medium increases.With the increase of the soot particle size,the soot particle scattering ability is enhanced,and the intensity transmissivity through the soot aerosol medium reduces.While the dust particle size is larger(such as 1 μm),the soot particle optical properties change gently.

soot aerosol;mie scattering;monte carlo;forward scattering

10.16185/j.jxatu.edu.cn.2016.09.007

2015-12-23

總裝備部重點實驗室研究基金項目(9140C360302120C36136)

侯宏錄(1960-)，男，西安工業(yè)大學教授，主要研究方向為光電測試、信息融合、復雜系統(tǒng)工程建模及仿真與作戰(zhàn)效能評估.E-mail：hlhou@sina.com.

O436.2

A

1673-9965(2016)09-0719-07