激光在海洋型气溶胶中传输的蒙特卡罗模拟

2018-12-27武文远龚艳春杨云涛

物理与工程 2018年6期

冷坤章曦武文远龚艳春杨云涛

(解放军陆军工程大学基础部，江苏南京 211101)

1 蒙特卡罗方法的提出

随着激光武器、激光雷达和激光通信等领域的发展，激光技术在国家安全防御方面发挥着越来越重要的作用，激光在大气中的传输特性一直是研究的热点，而激光在海洋大气中的传输特性研究相对较少。当激光在海洋大气中传输时，会受到海洋型气溶胶粒子、云雾粒子和雨滴等的散射和吸收，使激光能量在传输过程中不断衰减，影响激光系统的作用效能。其中，海洋型气溶胶粒子对激光传输的影响尤为显著，需要引起足够的关注。

1983年国际气象和大气物理协会提出的无云大气气溶胶模型将气溶胶分为6种：(1)水溶型粒子；(2)沙尘型粒子；(3)海洋型粒子；(4)煤烟；(5)火山灰；(6)75%硫酸可溶液滴。其中前4种出现在对流层，后两种出现在平流层[1]。这些粒子半径范围约为0.001～25μm，粒子数密度约为10～107cm-3。其中，海洋型气溶胶粒子的组成成分主要为水溶型粒子以及聚集模态和粗模态的海盐粒子[2]。

激光在海洋型气溶胶中传输特性的问题，实质上是一个存在多重散射的介质的辐射传输问题。由于激光在海洋型气溶胶中的辐射传输模型非常复杂，很难得到完全的解析式；实验手段又耗费巨大，也很难模拟所有情况。因此，蒙特卡罗方法作为一种以概率模型为基础，利用随机数和概率统计来进行研究的数值方法，是较为实用的一种手段。蒙特卡罗模拟光子输运过程有3种方法：直接模拟法、权重法、统计估计法。直接模拟法模拟过程简单清晰，但精度较差；权重法与统计估计法通过在模型中加入权重因子来减少计算误差[3]。而王红霞等人正是基于统计估计法研究了能见度、传输距离以及不同类型的气溶胶对激光传输透过率的影响[4]。其中，激光器自身的影响因素，光束初始半径、发散角等并没有考虑，并且光子到达接收平面处的空间分布也无法得到。本文基于权重法，结合光子出射时的高斯抽样和限制光子发散范围，提出一种新的蒙特卡罗模型，研究了传播距离、能见度和激光束初始均方根半径对到达接收平面处光子空间分布的影响。

2 理论分析

2.1 海洋型气溶胶的激光衰减特性分析

通常，不同地区气溶胶粒子的类型和浓度都不一样。一般情况下，各种组分的气溶胶粒子的谱分布可以用正态对数分布谱较好地描述[2]：

(1)

式中，Ni是单位体积内组分i的粒子总数；两个特征参数rmodN,i和σi为平均半径和标准差。

计算气溶胶的衰减系数时，通常将气溶胶粒子假设为球形粒子，再通过Mie散射理论[2]计算得到。但是实际的气溶胶粒子很多都不是球形，所以在计算非球形粒子的衰减系数时，可以采取T矩阵或者离散偶极子近似法(DDA)的方法[5-7]。本文假设海洋型气溶胶粒子为球形粒子。基于Mie散射理论计算气溶胶粒子的衰减系数时，还要利用到不同波长下气溶胶的复折射率参数。表1中给出了海洋型气溶胶粒子尺度谱的参数[4]以及在不同波长下的复折射率[8]。针对同一波长，空气相对湿度不同时，海洋型气溶胶复折射率变化较大，本文采用的是相对湿度为50%情况下的复折射率参数。

表1 海洋型气溶胶尺度分布模型参数及复折射率

计算气溶胶对激光传输的衰减系数，可以利用衰减截面、吸收和散射截面的方式求出，也可以利用能见度与衰减系数的关系近似求出[9]。对于尺度分布一定的球形气溶胶粒子，其衰减，散射以及吸收系数为[4]

(2)

式中，N为气溶胶粒子数密度；n(r)是气溶胶粒子尺度分布概率密度函数；σt，σs和σa分别为单个粒子的衰减、散射和吸收截面，其值可以根据Mie散射理论求出。根据Mie理论，衰减、散射、吸收截面如下定义[2]：

其中，an和bn为Mie系数；k为波数2π/λ。粒子的不对称因子表示前后散射的不对称程度，计算公式如下[1，4]：

其中，Qs为散射效率因子(Qs=σs/πr2)。单次散射率为ω=μs/μt。

气溶胶的数密度与能见度之间可通过关系式联立。对于人眼最敏感的波长(0.55μm)，能见度与衰减系数的关系如下：

(8)

根据上式可知，气溶胶粒子数密度为

(9)

把N代入公式(2)中，则可得能见度与1.06 μm波长下气溶胶的衰减系数之间关系为

(10)

表2为根据上述各式，在能见度为V=1km的情况下，计算得到的海洋型气溶胶粒子的消光参量，并与文献[6]中的参数进行了对比。

表2 海洋型气溶胶消光参量计算结果比较

2.2 蒙特卡罗计算模型

2.2.1 初始光束

根据激光谐振腔理论，激光输出应为高斯光束。在垂直于激光束传输方向的平面上，其光强分布为

(11)

式中,γs为激光束的均方根半径;用来描述光束的径向分布;x和y分别为垂直于激光束传输方向平面上的横坐标和纵坐标。在蒙特卡罗方法的源抽样处引入了高斯分布抽样模型，使发出光子的统计特性满足高斯光束的特性[10]。

2.2.2 光子在海洋型气溶胶中传输

海洋型气溶胶的衰减系数μt从上面可知，光子每一步的几何路径长度为

l=-ln(εr1)/μt

(12)

式中εr1为[0,1]之间均匀分布的随机数。之后判断光子是否被吸收，取一个[0,1]之间均匀分布的随机数εr2，如果εr2<ω，则光子被散射，假设光子的初始权值W=1，发生散射后，光子的权值变为W=W·ω；反之，光子被吸收，终止对该光子的跟踪[11,12]。

光子散射后，新的传播方向由散射相函数决定，可用Henyey-Greenstein函数来近似表示其散射相位概率函数，有时为得到更高精度的散射相位概率函数，也可以采用分段加权采样的方法求得相函数[13]。本文选取Henyey-Greenstein函数，其形式如下[14,15]：

(13)

式中,g为不对称因子;θ为散射角。假设方位角φ在(0,2π)内均匀分布，可表示为

φ=2πεr3

(14)

(15)

在光子运动过程中，如果光子跑到了激光发散角β之外的范围，也停止对其的跟踪。如图1所示：

图1 激光在气溶胶介质中传输示意图

图1中，L为激光传输的距离，假设传输过程中气溶胶粒子分布均匀。直到光子到达接收平面或者光子的权值小于阈值(本文取10-6)，就停止对其运动历程的跟踪。

3 模拟计算结果与分析

如图1所示，在激光传输一段距离后，垂直于传播方向设置一个接收平面，用于研究海洋型气溶胶对高斯光束传输光场空间分布的影响。模拟出射光子总数为106个，发散角β为10-3rad，衰减系数μt、单次散射率ω、不对称因子g由表2给出，其中ω和g为不随能见度变化的量，μt与能见度呈反比关系。首先我们考虑传播距离对接收平面光子空间分布的影响，蒙特卡罗模拟结果如图2所示，其中能见度为V=10km、初始光束均方根半径γs=0.5m；之后我们考虑能见度对接收平面光子空间分布的影响，蒙特卡罗模拟结果如图3所示，其中传播距离L=1000m、初始光束均方根半径γs=0.5m；最后我们考虑初始光束均方根半径对接收平面光子空间分布的影响，蒙特卡罗模拟结果如图4所示，其中传播距离L=1000m、能见度为V=10km。上述模拟结果均为相对湿度为50%情况下计算得到。

图2 接收平面光子的空间分布传播距离分别为(a) 100m； (b) 500m； (c) 1000m； (d) 2000m

图3 接收平面光子的空间分布能见度分别为(a) 5km； (b) 10km； (c) 20km； (d) 50km

由图2可知，随着传播距离的增加，接收平面处的光子空间分布仍未偏离高斯分布，但是到达接收平面的光子数明显减少，这是由于传播距离越远，光子被散射的概率越大，导致接收到的光子数减少。由图3可知，能见度较小时，光子空间分布仍未偏离高斯分布，但到达接收平面的光子数较少；而能见度较大时，光子空间分布仍未偏离高斯分布，而到达接收平面的光子数较多，这是因为能见度较大时，气溶胶粒子数密度较少，光子被散射的概率较小，导致接收到的光子数较多。由图4可知，激光束初始均方根半径越大，光子空间分布越分散，但仍未偏离高斯分布。

为了进一步分析光子空间分布与各个变量间的定量关系，图5给出了到达接收平面光子数与传播距离、能见度以及激光束初始均方根半径之间的对应关系。图6给出了光子的相对能量径向分布，即以接收平面中心为圆心，统计不同半径R范围内的光子能量E，并进行归一化。这里的光子能量为接收光子总能量。

由图5可知，随着传播距离的增加，到达接收平面光子数呈近似线性减少的趋势；而随着能见度的增加，到达接收平面光子数一开始迅速增加，然后增加趋势变缓；这一结论与文献[6]结果吻合。随着激光束初始均方根半径的增加，到达接收平面光子数减少，但变化并不明显。

由图6(a)可知，对于同一传播距离，随着半径R的增加，相对径向能量一开始迅速增加，然后趋于稳定，说明能量分布并未偏离高斯分布；对于同一半径R，随着传播距离的增加，相对径向能量减小，这些都与图2相吻合。由图6(b)可知，对于同一能见度，随着半径R的增加，相对径向能量一开始迅速增加，然后趋于稳定，说明能量分布并未偏离高斯分布；对于同一半径R，随着能见度的增加，相对径向能量增大，这些都与图3相吻合。由图6(c)可知，激光束的初始均方根半径越大，所需要的接收装置尺寸越大。激光束初始均方根半径为100mm时，所需接收装置尺寸半径最少为180mm；激光束初始均方根半径为300mm时，所需接收装置尺寸半径最少为830mm；激光束初始均方根半径为500mm时，所需接收装置尺寸半径最少为2320mm；激光束初始均方根半径为700mm时，所需接收装置尺寸半径最少为4750mm。这里接收装置尺寸半径为接收光子能量达到99%时半径R的值。