基于FDTD算法仿真光纤导波模式研究

1、基于FDTD算法仿真光纤导波模式的研究第1章 绪 论1.1 研究背景及意义 自1873年麦克斯韦(Maxwell)建立电磁场基本方程以来，电磁波理论和应用的发展已经有一百多年的历史。目前，电磁波的研究已深入到各个领域，应用十分广泛，例如无线电波传播、光纤通信和移动通信、雷达技术、微波、天线、电磁成像、地下电磁探测、电磁兼容，等等。电磁波在实际环境中的传播过程十分复杂：例如各种复杂目标的散射，复杂结构天线的辐射，在波导和微带结构中的传播，实际通信中城市环境、复杂地形及海面对电磁波传播的影响，等等。具体实际地研究电磁波的特性有着十分重要的意义。实验和理论分析计算是相辅相成的重要手段。分析计算途径需

2、要结合实际环境电磁参数求解麦克斯韦边值问题，通常只有一些经典问题有解析解。应当说，解析解具有重要指导意义。然而，由于实际环境的复杂性，往往需要通过数值解得到具体环境下的电磁波特性。 1966年K.S.Yee首次提出了一种电磁场数值计算的新方法时域有限差分(Finite Difference Time Domain，FDTD)方法。对电磁场E、H分量在空间和时间上采取交替抽样的离散方式，每一个E(或H)场分量周围有四个H(或E)场分量环绕，应用这种离散方式将含有时间变量的麦克斯韦旋度方程转化为一组差分方程，并在时间轴上逐步推进地求解空间电磁场。Yee提出的这种抽样方式后来被称为Yee元胞。FDT

3、D方法是求解麦克斯韦微分方程的直接时域方法。在计算中将空间某一样本点的电场(或磁场)与周围格点的磁场(或电场)直接相关联，且介质参数己赋值给空间每一个元胞，因此这一方法可以处理复杂形状目标和非均匀介质物体的电磁散射、辐射等问题。同时，FDTD的随时间推进可以方便地给出电磁场的时间演化过程，在计算机上以伪彩色方式显示，这种电磁场可视化结果清楚地显示了电磁波传播的整个物理过程，便于分析和设计。1.2 FDTD的发展与应用经过四十多年的发展，FDTD已发展成为一种成熟的数值计算方法。在发展过程中，几乎都是围绕几个重要问题展开的，即数值稳定性、计算精度、数值色散、激励源技术以及开域电磁问题的吸收边界条

4、件等。数值稳定和计算精度对任何一种数值计算方法都是至关重要的。A.Taylor和M.E.Brodwin4利用本征值方法给出了直角坐标系下FDTD的空间步长与时间步长之间的关系。X.Min等5研究了存在边界条件时FDTD的稳定性问题。对于数值色散，与实际的物理色散不同，它是由电磁场量在空间和时间上的对波动方程作差分近似处理造成的。这种色散引起的误差造成在计算区域内传播的电磁波逐渐畸变67。K. L.Shlager 等8比较了二维和三维空间中几种正交网格算法的色散误差。当采用其他变形或非正交网格时，必须重新分析其数值稳定性和色散特性911，P.Monk 和 E.Suli12分析了不均匀长方体网格算

5、法的稳定性。激励源的设计和引入也是FDTD的一个重要任务。目前，应用最广泛的激励源引入技术是总场/散射场体系12。对于散射问题，通常在FDTD计算空间中引入连接边界，它将整个计算空间划分为内部的总场区和外部的散射场区，如图1-1。利用Huygens原理，可以在连接边界处引入入射场，使入射场的加入变得简单易行。图1-1开域电磁问题中，为了在有限的计算空间内模拟无限空间中的电磁问题，必须在计算空间的截断边界处设置吸收边界条件。吸收边界条件从开始简单的插值边界，已经发展了多种吸收边界条件。在早期得到广泛应用的是G.Mur13的一阶和二阶吸收边界条件，它是基于B.Engquist和A.Majda14的

6、单向波方程而提出的差分格式，在FDTD仿真区域外边界具有0.5%到5%的反射系数。目前应用最广泛的是J.P.Berenger15-17的分裂式完全匹配层，以及Z.S.Sacks等18和S.D.Gedney20的各向异性介质的完全匹配层，它们可使FDTD模拟的最大动态范围达到80dB。另一方面，为了更好的拟合研究对象的形状，克服台阶逼近带来的误差，D.E.Merewether19提出了柱坐标系下的网格剖分方法，R.Holland20提出了球坐标系下的网格剖分方法，P.Monk和E.Suli12提出了变网格步长方法，S.S.Zivanovic等21和P.Thoma等22提出了亚网格技术（即在一般区

7、域采用粗网格，在电磁场快变区域采用精细网格）。利用这些技术，可以更精确地模拟各种复杂的结构，适应各种复杂的介质，提高了复杂介质中数值计算的精度。时域模拟一般获得的是近场电磁信息，为了得到诸如天线方向图或散射体雷达散射截面之类的远场信息，必须获得计算区域以外的频域场或瞬态场。多位学者在这方面做了许多工作，发展了一种高效的时域近远场变换方法23-26。借助这种方法，可以实现由计算区域内近场数据到计算区域外远场数据的外推。目前，粗糙面散射的FDTD，传递函数在FDTD中的应用，周期介质、各向异性介质、色散介质和含有集中元件的FDTD，以及网络并行FDTD技术等方面也取得了很大进展。FDTD在迅速发展

8、的同时，也获得了非常广泛的应用。目前，它几乎被应用到了电磁场工程中的各个方面，例如：电磁散射、生物电磁计量学、辐射天线的分析、微波器件和导行波结构的研究、散射和雷达截面的计算、周期结构的分析、电子封装和电磁兼容的分析、核电磁脉冲传播和散射的分析、以及微光学元器件中光的传播和衍射特性的分析等。随着新技术的不断提出，其应用范围和成效正在迅速地扩大和提高。第2章 FDTD的基本原理Maxwell方程是描述宏观电磁现象的一组基本方程。这组方程即可以写成微分形式，又可以写成积分形式。FDTD方法由Maxwell旋度方程的微分形式出发，利用二阶精度的中心差分近似，直接将微分运算转换为差分运算，这样达到了在

9、一定体积内和一段时间上对连续电磁场数据的抽样压缩2.1 Maxwell方程和Yee元胞根据27中电磁场基本方程组的微分形式，若在无源空间，其空间中的媒质是各向同性、线性和均匀的，即媒质的参数不随时间变化且各向同性，则Maxwell旋度方程可写成： （2-1a） （2-1b）式中，E是电场强度，单位为伏/米（V/m）；H是磁场强度，单位为安/米（A/m）；表示介质介电系数，单位为法拉/米（F/m）； 表示磁导系数，单位为亨利/米（H/m）；表示介质电导率，单位为西门子/米（S/m）；表示导磁率，单位为欧姆/米（）。在直角坐标系中，（2-1）式可化为如下六个标量方程： （2-2） （2-3）这六个

10、偏微分方程是FDTD算法的基础。 K.S.Yee3在1966年建立了如图2-1所示的空间网格，这就是著名的Yee氏元胞网格。在FDTD中，空间上连续分布的电磁场物理量离散的空间排布如图2-1所示。由图可见，电场和磁场分量在空间交叉放置，使得在每个坐标平面上每个电场分量被磁场环绕，每个磁场分量也被电场环绕。这种电磁场的空间结构与电磁感应和电磁波传播的规律相符，在每一个网格单元都能满足法拉第感应定律和安培环流定律。各分量的空间相对位置也适合于Maxwell方程的差分计算，能够恰当地描述电磁场的传播特性。同时，电场和磁场在时间上交替抽样，抽样时间间隔相差半个时间步，使Maxwell旋度方程离散以后构

11、成显式差分方程，从而可以在时间上迭代求解，而不需要进行矩阵求逆运算。因此，由给定相应电磁问题的初始条件，FDTD就可以逐步推进地求得以后各个时刻空间电磁场的分布。图2-1Yee氏网格及其电磁场分量分布 引入如下的差分近似方法对（2-2）、(2-3)式中的六个偏微分方程进行了差分离散。令代表E或H在直角坐标系中某一分量，在时间和空间域中的离散可记为 （2-4）式中，、和分别是长方体网格沿x、y、z方向的空间步长，是时间步长，i、j、k分别是沿x、y、z方向的网格编号，n是时间步数。对关于时间和空间的一阶偏导数取中心差分近似，具有二阶精度，即 （2-5a） （2-5b） （2-5c） （2-5d）

12、2.2 FDTD的基本差分方程根据上述原则，可将（2-2）、（2-3）式离散为如下的差分方程形式： （2-6a） （2-6b） （2-6c） （2-6d） （2-6e） （2-6f）式中， （2-7a）， （2-7b）(2-6)式就是FDTD的基本差分方程组。2.3 数值稳定性2.3.1 时间离散间隔的稳定性要求 考虑时谐场情况 （2-3-1）这一稳态解是下面一阶微分方程的解： （2-3-2）用差分近似代替上式左端的一阶导数，上面方程变为 （2-3-3）定义数值增长因子q为 （2-3-4）代入（2-3-3)式得 （2-3-5）它的解为 （2-3-6）如果差分方程（2-3-3）式的解趋于解析解（

13、2-3-1）式，则，于是增长因子q应为 （2-3-7）上式表明数值稳定性要求在时间步，足够小时增长因子，则其充分条件为 （2-3-8）2.3.2 Courant稳定性条件从麦克斯韦方程可导出电磁场任意分量均满足齐次波动方程 （2-3-9）考虑平面波的解 （2-3-10）采用有限差分近似 （2-3-11）因此（2-3-6）式被离散为 （2-3-12）上式又可改写为 （2-3-13）其中用到了（2-3-6）式。亦即 （2-3-14）此式即为空间和时间的离散间隔之间应满足的关系，又称Courant稳定性条件。其中C是电磁波的相速度，若采用等间隔离散，即，则有。在具体的空间步长选择上，一般要满足数值色

14、散对空间离散间隔的要求： （2-3-15）其中是无色散介质中的波长。2.3.2 数值色散 FDTD方程组是对Maxwell旋度方程进行差分近似，在进行数值计算时，将会在计算网格中引起数字波模的色散，即在FDTD网格中，电磁波的相速与频率有关，电磁波的相速度随波长、传播方向及变量离散化的情况不同而改变。这种关系由非物理因素引起，且色散将导致非物理因素引起的脉冲波形畸变、人为的各向异性和虚假折射等现象。显然，色散与空间、时间的离散间隔有关，如下式所示：（2-3-16）式(2-11)是三维情况下在FDTD方法中的单色平面波数值色散关系的一般形式，它表明FDTD计算中波的传播速度与传播方向有关。式中、

15、分别是波矢量沿、方向的分量，是角频率，是被模拟的均匀介质中的光速。与数值色散关系相对应，在无耗介质中的单色平面波，色散解析关系是： （2-3-17）由式（2-3-16）可知，当式（2-3-16）中的、均趋于零时，它就趋于式（2-3-17）。也就是说数值色散是由于用近似差分替代连续微分而引起的，而且在理论上可以减小到任意程度，只要此时时间步长和空间步长都足够小，但这将大大增加所需的计算机存储空间和计算时间，并使累积误差增加。因此，在实际计算中要根据问题的性质和计算机的软硬件条件来选择合适的时间步长和空间步长。为获得理想的色散关系，问题空间分割应按照小于正常网格的原则进行。一般选取的最大空间步长为

16、，为所研究范围内电磁波的最小波长。由上分析说明，数值色散在用FDTD法分析电磁场传播中的影响是不可能避免的，但我们可以尽可能的减小数值色散的影响。2.4 PML边界吸收条件 完全匹配层1618（Perfectly Matched Layer, PML）是1994年由J.P.Berenger首先提出，并将其设置在FDTD计算区域截断边界处，用来吸收外向电磁波。Berenger假设将电磁场分量在PML介质中分裂，并分别对各个分裂的场分量赋以不同的损耗。这就相当于在FDTD区域截断边界外设置了一种特殊的非物理的吸收介质层，该层介质的波阻抗与相邻介质的波阻抗完全匹配，因而外向波将无反射地穿过分界面进入

17、PML。同时，由于PML为有耗介质，而且不依赖于外向波的入射角和波阻抗，即使为有限厚度，外向波在其中也会迅速衰减。在实际计算中，PML是目前一种很常用的吸收边界条件，有很好的吸收效果，其总的网格噪声能量是使用普通吸收边界条件时的1/107，可使FDTD模拟的最大动态范围达到80dB。在电磁场理论中，当一平面波从媒质参数为的介质l入射到媒质参数为)的介质2时，如果，y方向上两介质波阻抗相等时无反射波，波完全入射并且波进入介质2后将按的差值指数减小。设适当值就可以使y方向的波迅速衰减。同理，时，x方向的入射波无反射，入射波将按的差值指数衰减。在二维条件下，将计算区域分为如图所示的九个区域，其中FD

18、TD仿真区为真空媒质参数为（0,0,0,0）；其余八个部分为PML边界，如图2-1所示。图2-2PML边界参数设置波进入PML介质层后遇壁被反射回来并从新进入FDTD防真区，这时大小衰减为，其中a为M层的厚度。我们可以根据要求的衰减系数c（即上式中的e指数项)来确定参数，的大小，并且由于要满足PML与自由空间阻抗相等，则，由此，也确定了。为了避免参数突变带来的反射，计算中采用了渐变的参数设置。设置PML为n层反射系数Re，PML变化阶数为P，则，一，这里x为由内到外的层序数，在PML层区，要采用不同于仿真区的计算方程，将分为和 （2-4-1）并对式(2-3-18)取差分，可得到E、H的差分的计

19、算表达式。 2.5 激励源 FDTD仿真中另一个比较重要的问题是激励源的模拟，即选择合适的入射波形式并加入到迭代运算中进行FDTD的仿真分析。目前来看有常用的两种激励源，一种是时谐场源，一种是脉冲波源，下面分别加以介绍。（1）时谐场源 （2-5-1）上面即为一个时域的正弦谐波，开始时刻为0，要达到稳态传播，一般需要35个周期，式中为角频率。（2） 高斯脉冲源 （2-5-2）式子中是一个常量，标识了脉冲的宽度，当时，脉冲出现峰值。2.6 FDTD计算所需时间步的估计 为了使计算达到稳定,通常计算所需要时间步按照电磁波往返穿越FDTD计算区对角线35次来估计。若FDTD计算区总元胞数为，则对角线上

20、元胞约为 （三维)。按照Courant稳定条件，设计算中心区，即穿越对角线一次需要时间步为。总计算时间步约需步。对于二维情况则约为。或者说，计算时间步大约等于FDTD计算区对角线上元胞数目的1220倍。实际上，计算所需时间步还与散射体具体形状、结构有关。第3章 EastFDTD仿真介绍本文主要是利用EastFDTD仿真软件进行光纤导波模式的研究，根据光纤的结构和特性参数建立EastFDTD的仿真模型，记录仿真数据及仿真模式的图片。下面简单介绍仿真流程及细节设置。 使用EastFDTD进行数值仿真大致可以分为四个步骤： 建立工程、设置工程属性；建立仿真模型；建立光源；建立数据及数据后处理 ；3.

21、1 第一步：建立工程、设置工程属性运行EastFDTD，选择“文件”“新建”“EastFDTD工程文件”，设置EastFDTD工程名，选择工程文件保存的位置，指定文件保存类型，设置后点击确定。具体参数设置如图3-1所示。 图3-1-1新建工程文档选择菜单“模型”“属性”或在建模编辑窗口中的右键菜单中选择“属性”，弹出“基本参数设定”窗口如图3-2所示，对该工程文件的基本参数进行设置。 图3-1-2基本参数设置基本参数的设定要进行全局参数、边界参数、计算参数、计算模式和变量的设置，这些变量的设置要根据仿真模型的结构、光源的波长以及具体实际情况而定。这里不做详细的介绍，具体的设置请参考EastFD

22、TD的参数设置模块。3.2 第二步：建立仿真模型 建立仿真模型分为两个步骤：建立元件和建立材料。 选择菜单“模型”“新建” “新建元件”或单击工具栏上的“新建元件”图标，设置元件的参数。 建立一个直径为80nm的圆柱单元，具体参数设置如图3-2-1所示。 图3-2-1元件参数设置 选择菜单“模型”“新建” “新建材料”或单击工具栏上的“建材料”图标，设置材料参数。 图3-2-2材料参数设置设置完模型的两个属性后，开始建立模型结构。选择菜单“模型”“新建” “新建结构”或单击工具栏上的“新建结构”图标。选择前两步建立的元件和材料。可通过右侧的“材料属性”和“元件属性”对选定的材料和元件参数做相应

23、的编辑。 图3-2-3结构参数设置3.3 第三步：建立光源选择菜单“模型”“新建”“新建光源”或单击工具栏上的“新建光源”图标，设置光源的参数。光源的设置具体要根据仿真模型的实际情况。各种类型光源参数的意义请参见软件中的光源模块。 图3-2-4高斯脉冲参数设置3.4 第四步：建立数据记录及数据处理选择菜单“模型”“新建”“新建记录”或单击工具栏上的“新建记录”图标，设置记录的参数。各种类型记录参数的意义参见记录模块。本例计算透射系数和反射系数，可选用“截面场平均”记录器来记录透过的平均场和反射的平均场，Ey分量反射平均场记录具体参数设置如图3-2-5所示。 图3-2-5截面场平均记录器参数设置

24、所有参数设置完毕，确认后，选择“模型”“运算”“*线程运算”，开始进行FDTD求解计算。计算过程中可双击“计算任务管理窗口”中的任务查看实时场。计算结束后数据将被保存输出到第一步全局参数中指定的路径，默认保存在工程文件同目录下。利用数学工具及图形处理软件，对数据进行简单的处理，就可得到想要的数据处理结果。EastFDTD中其他参数的设置具体参见软件中的其他模块，本文中不作具体介绍。第4章 建立仿真模型及数据处理本文是利用EastFDTD软件仿真研究光纤的导波模式。其中光纤的结构参数为：线芯直径6,折射率为1.5，包层直径72，折射率为1.47。由光纤的结构参数可知光纤的传导模式为单模。仿真所用光源为中心波长为1550nm，半宽为50nm的高斯脉冲。由光纤的结构参数和光源的波长范围，根据数值稳定性条件则可具体设置仿真模型。本文中设置的仿真模型如图所示。4.1 工程的属性参数 根据光纤的尺寸以及所要设置的波段选择微米为长度单位，选择自动相关联，在默认情况下软件自动将时间单位与长度单位相关联。其他参数设置根据仿真实际要求设置如图4-1-1。图4-1-1全局参数设置 本例为光纤结构，在Z方向均为周期分布，所以在Z方向上选择周期边界，X、Y选择开放边界，并使用默认设置。具体参