大气湍流中光传播的数值模拟?
马保科1,2， 郭立新1 吴振森1
（1.西安电子科技大学，陕西西安 710071 2.西安工程大学，陕西西安 710048 ）
摘 要 光在大气湍流中传播时，受大气分子、气溶胶等粒子的相互作用，将发生光束扩展、漂移和相干性退化等大气湍流效应，这些因素严重影响了光波的远场特性。文章从大气湍流中光传播的理论研究入手，分析了如何构造较为合理的大气湍流相位屏。进而采用McGlamery算法，对Kolmogorov谱下的大气湍流随机相位屏进行了数值模拟，并分析了光波从发射机经湍流大气传播到达接收机时的远场变化特性。研究表明，大气湍流的存在对光的远场传播质量造成很大的影响，研究结果也为大气湍流中与光传播相关的工程应用及自适应光学技术的完善提供了参考。 关键词 大气湍流；McGlamery算法；相位屏模拟； 大气结构常数； 中图分类号 TP391 文献标识码 A
1 引言
大气湍流是一个相当复杂的随机媒质系统，虽然物理学界对湍流的研究已经历了相当漫长的历史，但因涉及的因素千头万绪，其间的相互作用和关系也错综复杂，人们对其物理本质至今未能做到较为清楚的认识。因此，光在大气湍流中传播问题的研究仍存在理论和实验上的挑战[1,2]。通常，当光在湍流大气中传播时，光束截面内包含着许多的大气漩涡，这些漩涡各自对照射到它的那一部分光束形成衍射作用，可导致光束的强度和相位随机变化，进而表现出光束扩展，大气闪烁和相位起伏等大气湍流效应，从而严重降低了接收机的接收效率。目前，突破大气湍流的影响仍是光在随机介质中传播所要解决的关键问题[3]。早在20世纪中期，苏联的Obukhov便采用Rytov平缓微扰法由实验反演湍流特征。在闪烁的饱和现象被发现之后，物理学界又将Markov近似引入求解光场的统计矩，研究大气湍流下的光场特性[1]。然而，在中等起伏条件下，目前仍没有找到很好的解析处理方法。由于数值模拟能够从光的传播过程出发，较为清楚地反映出所涉及问题的物理本质，因而成为研究湍流效应的主要方法[4]。本文采用McGlamery算法[5]，对Kolmogorov谱下的大气随机相位屏进行了数值模拟，进而结合Huygens-Fresnel原理，模拟了在有无大气湍流的情况下，接收机处光场的变化特性。
2 大气湍流中光的传播
在折射率为n的随机媒质中，一束波长?，波数为k（k?2??）的单色波的电场E由Maxwell波动方程来描述[1,4]
?2E?k2n2(r)E?2?(E??lnn)?0 (1) 收稿日期: 2010年3月14日 收到修改稿日期：2010年 月 日
基金项目：教育部科技重点项目（105164） 通信作者: 马保科（1972），男，在读博士，副教授，主要从事随机介质中波传播方面的研究。Email: baokema2006@126.com
1
r1?其中，n(r)?1?n1(r)，（n1()，且n1(r)?0）为湍流折射率，(1)式左端最后一项
反映了偏振特性，当波长远小于湍流的内尺度（??l0）时，此项可以忽略不计，有
?2E?k2n(r)2E?0 (2)
因而，电场E的任一分量波动方程为
?2E?k2n(r)2E?0 (3)
如果介质的非均匀尺度远大于波长，可认为只存在前向小角度散射而无后向散射。对如
图1沿z方向的传播进行旁轴近似处理，并将光场写为E?uexp(ikz)，得
?2u?u?u?2?2ik?k2(n2?1)u?0 (4)
?z?z22222其中，?????x???y为横向算符。(4)式还可进一步化简为
2?(??????ik?ikH)(?ik?ikH)u?ik?H,?u?0 (5) ?z?z??z?2其中，算子H?1?k?2?2，??(n?1)。如果折射率不随传播距离而变化（或变化较小）则(5)式中忽略交换算子项，得到如下的广义抛物型方程
?u?ik(H?1)u (6) ?z通常，算子H不同的近似可得不同的抛粅型方程，而算子H最简单的近似即为其Taylor展式
2H?1?（k?2???(n2?1)）/2 (7) 当湍流场的折射率满足n(r)?1?n1(r)时，(7)式化简为
H?1?2k?2?2??n1 (8)
将(8)式代入(6)式得
?ui??2?u?ikn1u (9) ?z2k或
?2u?2u?u??2ik?2k2n1u?0 (10) 22?x?y?z(9)式即为忽略后向散射而具有抛粅型近似的波动方程。如果仅考虑介质折射率起伏对场
的作用，则它的右边只需保留与折射率有关的第二项，这时其解反映了在光的传播方向上由于积分光学路径所导致的相位调制，即
zu(r,z)?u(r,z')exp(?ik?'n1(r,?)d??)???z? (11)
?u(r,z')exp(iS)(11)式中，如果折射率起伏引起的相位变化S足够小，则可将真空传播和介质相位调制看作是相互独立并同时完成的两个过程。这样，如图1可将连续随机介质分割为一系列厚度为?z的平行平板，位于平板前的光场根据(9)式的真空解传播至平板的后面，然后被该平板引起的相位调制；这个场再经同样的真空传播和相位调制传播至下一个平板，依次形成最终的光场。也即将光在湍流介质中的传播等效为光在真空中通过一系列
2
薄的相位屏后的传播。
图1 随机介质中光传播的相位屏模型
Fig .1 The phase screen model in the atmosphere turbulence
这样，通过第i个相位屏后的光场由(9)式得
?izi?u(r,zi)?exp???2dz??exp(iS(r,zi))u(r,zi?1) (12)
?2kzi?1?其中，S(r,zi)?k?n1(r,z)dz，由于S(r,zi)的随机特性，不可能利用解析的方法
zi?1zi获得上式的解，而只能利用数值方法。
3 大气湍流随机相位屏的数值模拟
湍流介质中光传播的特殊性就在于介质折射率是随机起伏的，数值模拟的一个关键
也在于如何构造合理的相位屏，来正确反映介质折射率的变化特性，通常所构造的相位屏必须满足[4,5]： （1）、相位屏所代表的平板厚度?z应足够小，以确保其对场的振幅没有明显影响而只影响相位，即
2 (13) ?z??/?n2其中，?n为湍流折射率起伏均方差。 （2）、光场的变化特性与相位屏的构造方法应无关，即平板厚度应大于湍流介质非均匀元尺度，即
?z?L0 (14)
这里，L0为湍流外尺度。
（3）要使光在平板内的传播满足几何光学近似，则Fresnel尺度lF应小于湍流内尺度l0，即
2?z?l0/? (15)
目前，大气湍流随机相位屏的模拟方法较多[6]。这里主要利用McGlamery算法来模拟大气湍流相位屏，算法步骤如下：
1）生成二维具有高斯分布的复随机数矩阵?k?。
3
2）根据大气折射率功率谱函数，生成二维的功率谱密度函数矩阵 。
3) 将功率谱密度函数矩阵求算术根再乘以复随机数矩阵?k?，得到一个相位均匀分布在
???,???，振幅受功率谱密度函数调制的复随机数矩阵，即
?Y(?)???K(?)???(?)? (16)
4）对矩阵?Y(?)?进行逆Fourier变换得到其空间域形式，即
?Y(l)??IFFT?Y(?)????Y(?)?exp(j2?l?)dl0L (17)
5）将矩阵?Y(l)?分解为实部和虚部。这样，每一部分均独立代表一种随机的大气湍流相位屏
?Screen1??Re?Y(l)?，?Screen2??Im?Y(l)? (18)
在模拟相位屏时，相关参数的选择为[4,5]：网格间距?x?l0/3或?x??L/?，屏间
2距?z?l0。 /?，网格长度N?x??L/?0或N?x?2D（D为发射机孔径尺寸）
610Structure Function Df(r) Simulant valueTheoretical value 10510410310-1100101 Distance r (m)
图2 模拟生成的大气湍流相位屏 图3 相位结构函数的理论与模拟值比较
(Fig.2 The simulant atmosphere turbulent phase screen, Fig.3 Compare the simulant and the theoretical value
of the phase structure function )
图2是利用McGlamery算法模拟生成的Kolmogorov谱下的大气湍流相位屏，由于通常单层相位屏上相位分布的重复性较高，这里我们采用了多屏叠加的思路来尽量减小相位的重复性[7]。模拟屏上的相位结构函数定义为
D?(r)???(??r)??(r)??2?[Screen(??r)?Screen(r)]r?1n2 (19)
n其中，同时，与KolmogorovScreen(r)表示相位屏上r处的相位，n为屏上的网格数。
大气湍流功率谱相对应的相位结构函数理论值为[2]
D?(r)?6.88(rr0)5/3 (20)
4
通常，大气湍流的统计特性可以用相位结构函数来描述，因此，可以将相位结构函数作为验证模拟的相位屏正确与否的判断标准。图3比较了相位结构函数的理论值与模拟值。其中，湍流的内、外尺度分别为l0?0.01m米和L0?1m米，大气相干长度r0?8mm。可见，在垂直于传播方向的平面内横向距离r?2米（即r?L05）的范围内，相位结构函数的理论值与Kolmogorov谱下的模拟值吻合较好，但当r?L05时，两者相差较大。因此，在湍流谱及上述参数选定的情况下，为了达到模拟大气湍流的准确性和避免大尺度湍流起伏所带来的误差，相位屏的尺寸最好选为L05?L05大小。
4 大气湍流下传播光场的数值模拟
图4 发射机平面上的模拟光斑
Fig.4 The stimulant facula on the transmitter plane
结合上述Kolmogorov谱下模拟生成的相位屏以及利用Huygens-Fresnel原理，下面就光在大气湍流中水平传播到达接收机处的光场进行数值模拟。这里，相关参数的选择为：发射光波为??0.6328?m的平面波,发射孔径为D?20cm的方形孔径，传播距离约10km, 湍流的内、外尺度分别为1cm和10m，传播路径上分布40个相位屏, 垂直于传播方向的平面被分成了100×100个小网格, 网格宽度选择为湍流的内尺度大小。
5