楊 穎，黃竹青，2，曹小玲

(1.長沙理工大學能源與動力工程學院，湖南長沙410110;2.湖南有色金屬職業(yè)技術(shù)學院，湖南株洲412006)

1 引言

隨著激光技術(shù)和計算機的飛速發(fā)展，光學法微粒粒度測量技術(shù)在越來越多的領(lǐng)域中出現(xiàn)，其中很大一部分是以光散射理論為基礎(chǔ)。Mie散射理論是研究各向同性的球形粒子光散射的基礎(chǔ)，又稱為大粒子彈性散射，其適應(yīng)范圍從幾納米到幾毫米甚至更大的氣體、液體和固體，如空氣中的一些氣體分子、硫氧化物分子、霧、云、汽輪機中的水滴等。針對球形粒子，本文利用Matlab軟件，選取合適Mie散射計算方法，對散射強度函數(shù)、相函數(shù)等進行了編程計算和分析。

2 球形粒子散射理論計算

Mie散射理論計算非常復(fù)雜，特別是Mie系數(shù)的計算，不少國內(nèi)外學者對算法進行了研究改進，并發(fā)表了若干種方法［1－5］，本文采用的是一種比較穩(wěn)定且精度高的算法，具體見2.2。

2.1 散射強度函數(shù)

根據(jù)Mie散射理論，當入射光是平面偏振光時，總散射光強的垂直分量Ir和平行分量Is表達式如下(φ是極化角):

總散射光強表達式為:

當入射光為非偏振光時總散射光強為:

其中，d是觀察點與微粒的距離;i1(θ)和i2(θ)是散射強度函數(shù);S1(θ)和S2(θ)是振幅函數(shù)，它們之間的關(guān)系如下:

由以上公式可知，總散射光強IS、散射強度函數(shù)i1(θ)、i2(θ)和振幅函數(shù) S1(θ)、S2(θ)都是 Mie系數(shù)an和bn、散射角分布函數(shù)πn和τn的無窮級數(shù)(關(guān)于截止階數(shù)nstop的確定在本節(jié)末給出)。所以只要計算出Mie系數(shù)an和bn、散射角分布函數(shù)πn和τn即可。

2.2 Mie系數(shù)的計算

Mie散射系數(shù)的計算方法中，用定義式計算比較困難，且非常耗時，所以一般采用基于定義式改進的向前或向后遞推算法。而當顆粒吸收性很大或者折射率虛部絕對值非常大時，向后遞推計算有可能超過計算機允許的范圍而溢出，因此本文采用向前遞推算法，根據(jù)文獻［6］對Mie系數(shù)的重新構(gòu)造進行計算，改進后的Mie系數(shù)表達式如下:

其中:

3)執(zhí)行情況檢查(C)：在成本控制執(zhí)行過程中，實時收集和掌握成本控制數(shù)據(jù)和信息，對成本控制執(zhí)行情況進行動態(tài)分析和前景預(yù)測，對執(zhí)行情況進行績效考核和激勵，動態(tài)和前瞻性地監(jiān)控成本管理的進展，確保項目投資在成本控制總目標內(nèi)進行。

其中，α=2πr/λ 是無因次粒徑參量;β=mα，m 是微粒相對于周圍環(huán)境的折射率。

求 Mie系數(shù) an和 bn轉(zhuǎn)化為計算 An(α)、Bn(α)、Ln(α)、Ln(β)。An(α)和 Bn(α)可以用向后遞推計算，遞推公式如下:

Ln(x)采用向前遞推法計算，遞推公式為:

由Lentz的連分式算法求初始值 Ln*(x)［7］，具體算法如下:

根據(jù)Lentz的證明:

通過計算比較Gnk與Gn(k－1)之差的絕對值滿足誤差要求，計算截止，再將Gnk代入公式(11)得到Ln*(x)。

散射角分布函數(shù)πn和τn同樣采用遞推方法求(其中θ是散射角):

初始值:π0=0;π1=1;π'0=π'1=0

圖1 截止階數(shù)的選取方法Fig.1 Selection method of cutoff number(r=2μm)

3 特性分析

散射光強的分布主要和無因次參量α、折射率m和散射角θ有關(guān)，散射(消光、吸收)系數(shù)和無因次參量α、折射率m有關(guān)，微粒群的散射特性還和分布模型有關(guān)，因此本文主要針對這幾個參數(shù)進行分析。而入射光波長λ的變化相當于無因次參量α的變化，所以改變波長相當于改變α。

3.1 散射強度函數(shù)

單個微粒的散射強度函數(shù)i1(θ)和i2(θ)表達式如公式(4)。結(jié)合Mie系數(shù)和角度分布函數(shù)進行編程計算，得到散射強度的分布關(guān)系，如圖2，前向散射明顯強于后向散射，且微粒半徑越大越明顯，散射光強大部分集中在前向小角度范圍內(nèi)，這被稱為米氏效應(yīng);圖2(a)為不同微粒大小之間的對比，散射光強呈現(xiàn)振蕩現(xiàn)象，隨著微粒半徑的增大，振蕩越明顯。圖2(b)為不同折射率之間的對比。由于耗散介質(zhì)的吸收，散射光強明顯比非耗散介質(zhì)弱，振蕩現(xiàn)象隨著折射率虛部絕對值的增大逐漸消失，變得平滑，尤其是后向散射。

圖2 散射強度隨角度分布關(guān)系Fig.2 The distribution of the scattering intensity varies with the Angle

3.2消光系數(shù)、散射系數(shù)和吸收系數(shù)

在非耗散介質(zhì)中，折射率是一個實數(shù)，表示沒有吸收;在耗散介質(zhì)中，復(fù)折射率m=n+iη，n反映了介質(zhì)的散射特性，η反映了介質(zhì)的吸收特性。消光系數(shù)和散射系數(shù)的計算公式分別如下:

通過Matlab對公式(14)進行計算得到各個系數(shù)和微粒大小的關(guān)系，如圖3。非耗散介質(zhì)比耗散介質(zhì)系數(shù)曲線振蕩明顯，曲線有一系列的極大值和極小值，并出現(xiàn)微小的毛刺現(xiàn)象，折射率實部越大，毛刺現(xiàn)象也越明顯，極大值和極小值在同一粒徑范圍內(nèi)增多，消光系數(shù)趨近于2的速度越快;而對于折射率虛部不為0的耗散介質(zhì)，系數(shù)隨粒徑變化的振蕩減弱，毛刺現(xiàn)象消失，曲線變得平滑。當折射率虛部增大，即微粒的吸收特性增強時，系數(shù)曲線變得非常光滑。圖3(c)中給出了消光系數(shù)，散射系數(shù)，吸收系數(shù)的關(guān)系，當顆粒為非耗散介質(zhì)時，吸收系數(shù)Qabs=0，Qext=Qsca;當顆粒為耗散介質(zhì)時，Qext=Qsca+Qabs。

圖3 各系數(shù)和微粒大小的關(guān)系Fig.3 Relation of the each coefficient and the particle size

3.3 散射相函數(shù)

光在某個給定方向單位立體角中散射的能量和所有方向上平均單位立體角中的散射能量之比，稱為散射相函數(shù)。它表示距離粒子很遠處的散射能量隨角度分布狀況［9－11］。當入射光為平面偏振光時，單個水滴相函數(shù)的表達式為:

對于單個微粒，散射系數(shù)相當于一個常數(shù)，其特性與散射強度函數(shù)相似。圖4給出了單個水滴的散射相函數(shù)和散射角的關(guān)系。

圖4 單個粒子的散射相函數(shù)Fig.4 Scattering phase function of single particle

在研究非單一分散系的散射特性時，假設(shè)微粒群中的各個微粒是相互獨立的散射體，即各個粒子之間的散射光互不相關(guān)，所以微粒群的散射特性是各個微粒散射特性的線性疊加。本文以適合描述破碎形成的微粒群的Rosin－Rammler分布函數(shù)表示微粒群尺寸分布，其尺寸分布概率密度函數(shù)表達式為:

其中，K是分布參數(shù);珋r是質(zhì)量中間半徑，它表示大于質(zhì)量中間半徑的顆粒數(shù)占50%，散射相函數(shù)與微粒群的數(shù)密度無關(guān)［12］，得到粒子群的相函數(shù)表達式為:

散射相函數(shù)滿足歸一化條件:∫P(θ)dθ=1

又有散射能量累積分布函數(shù):

表示散射角小于φ的散射能量占所有散射光的比例，因此有:F(π)=1。

圖5 微粒群散射相函數(shù)Fig.5 Scattering phase function of particle swarm

通過軟件對公式(17)進行數(shù)值積分，得到微粒群散射相函數(shù)與散射角的關(guān)系，如圖5。微粒群的光散射分布在前向比較平滑，后向小角度出現(xiàn)振蕩，與圖4單個球形微粒的散射特性比較，振蕩明顯減弱。分別對不同的分布參數(shù)進行計算、對比，不同的分布參數(shù)對散射相函數(shù)的影響比較小，如圖5(a);不同微粒大小對散射光的分布的影響如圖5(b)，與單個微粒相比，不同質(zhì)量中間半徑之間的差別減小，特別是在前向30°～60°之間。

4 結(jié)論

利用改進后的遞推方法求Mie系數(shù)，使計算速度大大的提高，而采用Lentz的連分式方法使計算結(jié)果更加精確。針對球形微粒，入射光為平面偏振光，本文分別對單個微粒和微粒群的散射特性進行了計算分析，對比各個參數(shù)對散射強度函數(shù)、散射相函數(shù)、散射(吸收、消光)系數(shù)等的影響，對光散射實驗研究具有指導作用，有助于對實驗數(shù)據(jù)的分析。

［1］ V DAVE.Scattering of visible light by large water sphere［J］.Appl Opt，1969，8(2):155 －161.

［2］ Wang Shimin，Zhu Zhen，Ye Mao，et al.Two improved numerical algorithms of mie theory in the particle sizing by light scattering method［J］.Acta Metrologica Sinica，1999，20(4):279 －285.(in Chinese)王式民，朱震，葉茂，等.光散射粒度測量中Mie理論兩種改進的數(shù)值計算方法［J］.計量學報，1999，20(4):279－285.

［3］ Wang Lianfen，Liu Guangbin，Zheng Zhong，et al.Onwards continued fraction algorithm for Mie scattering and numerical simulation［J］.Optical Technique，2009，35(4):541 －543.(in Chinese)王蓮芬，劉光斌，鄭重，等.Mie散射向前遞推連分式算法及其數(shù)值模擬［J］.光學技術(shù)，2009，35(4):541－543.

［4］ Wang Shaoqing，Ren Zhongjing，Zhang，Ximing，et al.New algorithm of mie scattering parameters［J］.Laser Journal，1997，18(3):9 －13.(in Chinese)王少清，任中京，張希明，等.Mie散射系數(shù)的新算法［J］.激光雜志，1997，18(3):15 －22.

［5］ 蔡小舒，蘇明旭，沈建琪，等.顆粒粒度測量技術(shù)及應(yīng)用［M］.北京:化學工業(yè)出版社，2010.

［6］ Shen Jianqi，Liu Lei.An improved algorithm of classical Mie scattering calculation［J］.China Powder Science and Technology，2005，11(4):1 －5.(in Chinese)沈建琪，劉蕾.經(jīng)典Mie散射的數(shù)值計算方法改進［J］.中國粉體技術(shù)，2005，11(4):1 －5.

［7］ W J Lentz.Generating bessel funtions in Mie scattering calculations using continued fractions［J］.Appl Opt，1976，15:668 －671.

［8］ Wiscombe W J.Improved Mie scattering algorithms［J］，Applied Optics，1980，19:1505 －1509.

［9］ Zhang Heyong，Ren Deming，Zhao Weijiang，et al.Research of scattering phase function of spheric particle under polarized state［J］.Opto － Electronic Engineering，2008，35(9):70 －75.(in Chinese)張合勇，任德明，趙衛(wèi)疆，等.偏振狀態(tài)下球形粒子的散射相位函數(shù)研究［J］.光電工程，2008，35(9):70－75.

［10］ Yuan Yijun，Ren Deming，Hu Xiaoyong.Computing scattering phase function by recursive formula of Mie Theory［J］.Chinese Journal of Light Scattering，2005，17(4):366 －371.(in Chinese)袁易君，任德明，胡孝勇.Mie理論遞推公式計算散射相位函數(shù)［J］.光散射學報，2005，17(4):366 －371.

［11］ Du Zhufeng，Lu yimin，Huang Tiexia.The approximate calculation of scattering phase function［J］.Journal of Optoelectronics·Laser，1998，9(5):422 － 425.(in Chinese)杜竹峰，盧益民，黃鐵俠，等.散射相位函數(shù)的近似計算［J］.光電子·激光，1998，9(5):422 －425.

［12］ Huang Zhuqing.Investigation of exhaust wetness measurement method in steam turbine based on mie scattering of laser［D］.Changsha:National University of Defence Technology.2010.1.(in Chinese)黃竹青.基于光散射的汽輪機濕蒸汽特性測量方法研究［D］.長沙:國防科技大學，2010.