基于光线光学的非线性自聚焦现象的仿真分析*

阮望超 岑兆丰 李晓彤 刘洋舟 庞武斌

(浙江大学现代光学仪器国家重点实验室,杭州 310027)

(2012年5月7日收到;2012年9月21日收到修改稿)

1 引言

对非线性自聚焦的研究始于1962年[1],至今已有许多研究自聚焦的理论[2−4],并在一定程度上得到了实验[5−7]的验证.自聚焦的仿真方法可分为波动光学方法和光线光学方法.前者是数值求解非线性薛定谔方程,该方法又可以分为两类[8],即有限差分法和伪频谱法.目前被广泛用于解非线性色散介质的脉冲传输问题的方法是分步傅里叶方法(split-step Fourier method,SSFM)或称光束传输法(beam propagation method,BPM)[9],它属于一种伪频谱法.此外,许多BPM的改进形式已经被提出来[10,11],提高了计算速度等.早期的数值研究是基于近轴近似的标量波理论的,因此得出光束崩塌的结论.后来人们引入了非近轴修正并考虑到光束的矢量特性[12−14],证明了光束并没有崩塌,而是继续传播并形成多个焦点.在国内,中国科学院上海光机所的范滇元、文双春等[15,16]对非线性自聚焦效应进行了深入的研究,利用BPM方法编制了仿真程序并发表了许多相关论文.波动光学方法可以同时考虑自聚焦、群速度色散、衍射、损耗[17]等效应,因此求解过程相对复杂,通常可以采用一些近似来简化问题,例如近轴近似、假设光波的传播方向一直沿某一坐标轴方向(如z方向)、解析结果事先假定解的形式(隐含了不同位置的光分布自相似的前提).光线光学方法最早由Censor[18]提出来,他设定光强与光线间距成反比,当光线相交时就会得到光强极大的结果.Kasparian等[19]提出追迹大量光线的方法来仿真自聚焦效应,并认为在一个相当小的步长内折射率分布是均匀的,当激光光强较大时自聚焦点位置的计算结果比较准确,但他没有涉及非近轴区域光的传输和对网格统计的光强进行减少量化误差的处理.对于光线光学方法中的光强,除了可以使用网格来统计,还可以用Wigner分布函数(Wigner distribution function,WDF)[20]来求得,但是非近轴的情况下求解WDF的传输方程比较困难.四川大学张玉成和吕百达[21]推导了自由空间中部分相干非近轴光的WDF传输方程,但是没有用于非线性自聚焦.Gao等[22]利用WDF和迭代光追方法仿真自聚焦过程,可以用光线直观地表示激光传输的过程,但没有提及是否采用近轴近似.中国工程物理研究院激光聚变研究中心引入光强功率谱密度表示归一化光强随频率的分布,并利用“星光II”激光装置观察到了激光在硅酸盐介质中传输的非线性自聚焦现象[23].

为了对含有非线性自聚焦介质的实际系统进行仿真,本文提出了一种改进的光线追迹方法.在不引入近轴近似且考虑波像差的前提下,在光传输路径上垂直于光轴进行切片,在每一切片上建立正方网格,并对网格统计的光强进行减少量化误差的恢复处理,利用其结果计算光强对折射率的调制作用,而相邻两切片间只考虑光的线性传输,其结果能够在宏观上直观地体现出激光的自聚焦传输过程.

2 原理与算法

2.1 非线性光线追迹

自聚焦是强激光在非线性介质中传输时的一种非线性现象.在非均匀线性介质中,假定在很小的步长下光所走过的区域介质是均匀的,从而可以用数值求解光线方程的方法来仿真光传输过程;在非线性自聚焦情况下,还需要考虑光强对介质折射率的调制作用.本文考虑的光强对折射率的作用是瞬时响应的强局域作用.

在线性的非均匀介质中,光线追迹可以通过求解光线方程

来达到,其中n为折射率,s为光程.这本质上是求解二阶偏微分方程,典型的解法有泰勒展开法、欧拉法、龙格库塔法等[24].一般光学软件的光线追迹常用龙格库塔法.在非线性自聚焦介质中,由于必须考虑光强对折射率的影响,则上述方法会有一些问题:龙格库塔法属于高阶单步法,在计算某一位置的参数时,需要计算几个辅助的位置处的参数,这些位置是在当前位置之后的,光还没有追迹到那些位置上,在非线性情况下这些辅助的参数没办法直接求得,故该方法不能直接用于非线性自聚焦的光线追迹.

本文将显式亚当斯法用于求解光线方程,不会产生上述问题,且不需要迭代,而隐式的亚当斯法或者修正的亚当斯法则需要迭代,计算量较大.

光学方向余弦为

定义:

则(1)式可以改写为

给定初始值R0,T0及n的分布就可以进行逐步计算.

对于初值问题

显式亚当斯法的计算式为

其中上标的一撇表示对x求导.

由于只有一个初值,需要用多步法(如龙格库塔法)先算出几个初值(不同阶数的算法需要的初值数量不一样),可以使用迭代的方法提高初值的精度或者使用较小的步长追迹得出精度较高的初值.

光追过程中采用正方形网格来统计光强,并依此算出折射率及折射率梯度分布,最后算出折射率梯度对光线偏折的影响,从而得到光线的进一步传输,如此循环以至像面.光线追迹程序算法流程如图1.

在程序中,由于需要将所有光线都追迹到同一个面再进行光强统计,所以步长h必须根据切面距离∆z来自适应,即

图1 光线追迹主程序流程图

2.2 光强分布量化误差的抑制

用统计光线的方法求光强时,由于网格有一定大小和边界,因此必然引入量化误差.这种误差的特点是均值为零,且其空间频率是可估计的,一般比光强分布的空间频率高.

图2 抑制量化误差前后的光分布

一般的光学设计、仿真软件在处理量化误差时有两个方法,一是增加光线数和网格数,以减少量化误差;二是对光强分布进行平滑,采用高斯平滑算子进行多次平滑.但是这两个方法都有各自的缺点,前者会增加计算机的时间复杂度和空间复杂度,而且也不能完全消除量化误差;后者会导致过度平滑(收敛于高斯分布).

本文采用一种基于梯度信息的算法,对含有量化误差的光分布进行处理,得到相对光滑的光分布,见图2.

3 数值结果分析

强激光在大气中传输时,会由于非线性自聚焦而形成一些细丝.折射率受光强的调制关系是n=n0+n2I,其中I为光强,n2是非线性系数,n0是不考虑非线性情况下的折射率.自聚焦的临界光功率为实际光功率P＞P时出现自聚cr焦现象.空气中 n2=3×10−23m2/W,λ=800 nm,Pcr≈3 GW.

在程序实现上,采用m×n的网格,为了便于表述,令m=200,n=1,此时相当于一维的光分布.每一个网格内有20—30条光线,因此对于光束宽度在厘米量级的激光,初始时刻的横向分辨率约5µm.步长∆z在毫米量级.

3.1 轴上光强变化

现有的理论和实验已经证明,强激光在传输时会经历多次自聚焦的过程,在不同的传输距离上会形成多个自聚焦点.在仿真程序中,我们得到轴上光强随传输距离的变化呈现出近周期性的变化(见图3).

图3 轴上光强随传输距离的变化(P=100P cr)

3.2 整体光束的传输过程

光束的传输如图4所示,横坐标表示光束横截面半径(−15—15 mm),纵轴表示传输距离(0—100 m).不同横截面位置由于其光强梯度不一样(相应的折射率梯度也不一样),光线会聚的程度不一样,梯度大的位置光线会聚更快,因此形成了一股股的局部光束(二维情况下表现为一个个半径不等的圆环).这种情况类似线性介质中透镜的像差,文献[25]也曾指出在第一个自聚焦点到达之前,不同孔径带的光会聚点不一样会造成多个自聚焦焦点这种类似像差的现象.

图4 光线轨迹图

在另一个实验条件下,这种分散的细光束表现得更明显,如图5所示.在二维情况下,相应的表现是形成一系列的环,实验上可以观察到这种现象[26].

图5 分散的细光束(图中箭头所指)

在采用近轴近似的仿真程序中,常常可以见到类似梯度折射率介质中的重复自聚焦的现象,如图6(a).在我们的程序中,选取靠近光轴部分的光,如图6(b),则可见在近轴部分光传输也是类似梯度折射率介质中的重复自聚焦过程.

图6 (a)梯度折射率介质中的光传输;(b)非线性自聚焦近轴光线

由于本文中的方法追迹的是实际光线,不引入近轴近似,因此可以得到近轴区以外区域自聚焦及成丝(环)的情况,这对于强激光系统安全是有重要意义的.

3.3 第一个焦点的位置

对于不同光功率的激光,其自聚焦的起始位置是不一样的.我们对各种功率的激光进行仿真,得出了第一个自聚焦焦点位置随光功率的变化曲线,如图7.与Wanger等[27]基于波动光学方法推导出的结果进行对比,可见我们的结果在光功率较大时是足够准确的.

图7 不同光功率下第一个自聚焦焦点的位置

4 结论

采用几何光学的方法,通过光线追迹来仿真非线性自聚焦现象,有效避免了波动光学方法中的近轴假设、自相似等近似,从而在宏观上直观地体现实际光传输的过程.通过采用在光传输路径上垂直于光轴切片的方法把非线性光传输转化为切片端面上光强对折射率的非线性调制和切片间的线性光传输;切片间的传输采用了亚当斯法求解光线方程,很好地解决了龙格库塔法无法用于非线性光线追迹的问题;此外,还对网格统计的光强的量化误差进行了有效的抑制.对轴上光强随距离变化的仿真结果表明,强激光的非线性自聚焦过程会产生多个焦点,而第一个焦点的位置随着入射光功率的增大而更加靠近入射位置,准确地获得第一个焦点的位置对于强激光远距离传输有重要的意义;从整体光束的传输过程来看,光束横截面上不同位置的光强梯度(折射率梯度)不一样会导致圆环状光分布(一维情况下是传输方向偏离整体光束传输方向的局部细光束).基于光线追迹的非线性自聚焦仿真方法与其他方法相比有一个明显的优点,即可以把非线性介质和线性介质结合起来进行光传输的仿真.目前基于光线追迹的光学设计和仿真软件的分析功能已很齐全,若在其中添加非线性自聚焦的传输模块,即可仿真分析实际强激光传输系统.

[1]Askar’yan G A 1962 Soviet.Phys.JETP 15 1088

[2]Chiao R,Garmire E,Townes C H 1964 Phys.Rev.Lett.13 479

增量式光电编码器是将一系列脉冲通过附加电路处理得到数字量的编码器。光电编码器的码盘一般多用玻璃材料制成，玻璃表面镀上一层不透光的金属薄膜，然后在薄膜上刻制码道，数量从几百条到几千条不等。这样整个码盘圆周被等分成n个透光的槽。

[3]Kelley P L 1965 Phys.Rev.Lett.15 1005

[4]Shen Y R(Translated by Gu S J)1987 The Principles of Nonlinear Optics(1st Ed.)(Beijing:Science Press)p327(in Chinese)[沈元壤著(顾世杰译)1987非线性光学原理(上册)(北京:科学出版社)第327页]

[5]Hercher M 1964 J.Opt.Soc.Am.54 563

[6]Braun A,Korn G,Liu X,Du D,Squier J,Mourou G 1995 Opt.Lett.20 73

[7]Liu W,Kosareva O,Golubtsov IS,Iwasaki A,Becker A,Kandidov V P,Chin S L 2003 Appl.Phys.B 76 215

[8]Agrawal G P(Translated by Jia D F,Yu Z H)2002 Nonlinear Fiber Optics,Third Edition and Applications of Nonlinear Fiber Optics(1st Ed.)(Beijing:Publishing House of Electronics Industry)p33(in Chinese)[阿戈沃著(贾东方,余震虹译)2002非线性光纤光学原理及应用(第一版)(北京:电子工业出版社)第33页]

[9]Hardin R H,Tappen E D 1975 SIAM Rev.Chronicle.15 423

[10]Sinkin O V,Holzl¨oner R,Zweck J,Menyuk C R 2003 J.Lightwave Technol.21 61

[11]Godoy-Rubio R,Romero-Garc´ıa S,Ortega-Mo˜nux A,Wang¨umert-P´erez J G 2011 J.Opt.Soc.Am.B 28 2142

[12]Feit M,Fleck J 1988 J.Opt.Soc.Am.B 5 633

[13]Chi S,Guo Q 1995 Opt.Lett.20 1598

[14]Li Z D,Guo Q,Lin WG 1999 Chin.J.Lasers 26 711(in Chinese)[李忠东,郭旗,林为干1999中国激光26 711]

[15]Wen S H 2001 Ph.D.Dissertation(Shanghai:Graduate University of Chinese Academy of Sciences)(in Chinese)[文双春2001博士学位论文(上海:中国科学院研究生院)]

[16]Wen S C,Fan D Y 2000 Acta Phys.Sin.49 460(in Chinese)[文双春,范滇元2000物理学报49 460]

[17]Wen S C,Fan D Y 2000 Acta Phys.Sin.49 1282(in Chinese)[文双春,范滇元2000物理学报49 1282]

[18]Censor D 1977 Phys.Rev.A 16 1673

[19]Kasparian J,Solle J,Richard M,Wolf J P 2004 Appl.Phys.B 79 947

[20]Dragoman D 1996 Appl.Opt.35 4142

[21]Zhang Y C,L¨u B D 2004 Opt.Lett.29 2710

[22]Gao H H,Tian L,Zhang B,Barbastathis G 2010 Opt.Lett.35 4148

[23]Peng Z T,Jing F,Liu L Q,Zhu Q H,Chen B,Zhang K,Liu H,Zhang Q Q,Cheng X F,Jiang D B,Liu H J,Peng H S 2003 Acta Phys.Sin.52 87(in Chinese)[彭志涛,景峰,刘兰琴,朱启华,陈波,张昆,刘华,张清泉,程晓峰,蒋东镔,刘红婕,彭翰生2003物理学报52 87]

[24]Qiao Y T 1991 Gradient Index Optics(1st Ed.)(Beijing:Science Press)p170(in Chinese)[乔亚天1991梯度折射率光学(第一版)(北京:科学出版社)第170页]

[25]Dyshko A L,Lugovoi V N,Prokhorov A M 1967 ZhETF Pis’ma 6 655

[26]Garmire E,Chiao R Y,Townes C H 1966 Phys.Rev.Lett.16 347

[27]Wanger W G,Haus H A,Marburger J H 1968 Phys.Rev.175 256