第四章光束传播法基础第一节数值计算方法1．电磁场数值计算它是一种基于麦克斯韦方程组，建立逼近实际工程电磁场问题的连续型的数学模型，（合理的假设）然后采用相应的数值计算方法，经离散化处理，(合适的方法，使离散化的模型既能反映连续型模型的特性，又便于计算机分析)把连续型数学模型转化为等价的离散型数学模型，计算出待求离散数学模型的离散解(数值解)，从而获得相应结果的一种方法。
2.数值方法分类:时域分析、频域分析。
时域分析：模拟光在波导中的传播过程频域分析：求解波导模式时域分析逼真：把原来因为速度太快、结构太小、不可见的现象模拟出来，能够直观地展示。
求解：波导连接、耦合、非线性特性、波导模式。
频域分析：光场分布、给定具体结构波导的模式的有效折射率(色散、偏振)、损耗(材料吸收、结构本身导致)等。
问题: 频域结果能否推得时域信息? 反之?3．常用数值方法简介(1)有限差分法(频域有限差分法)（20世纪50年代出现）利用划分网格的方法将定解区域离散化为网格离散节点的集合，然后基于差分原理，以各离散点上函数的差商来近似替代该点上的偏导数，这样待求的偏微分方程定解问题可转化为一组相应的差分方程的问题。
根据差分方程组，解出各离散点上的待求函数值，即为所求定解问题的离散解，再应用插值方法便可从离散解得到定解问题在整个场域上的近似解。
原理:偏导→差分方法特点：原理简单、通用性好；对复杂结构，计算量大(矩阵运算)。
(频域分析)适用范围：计算光波导的模式求解。
现状：适用于较简单结构的分析。
但有限差分(偏导→差分)法广泛应用于数值方法中(2) 有限元法20世纪40年代提出,其在电磁问题方面的应用有约40多年历史。
以变分原理为基础，把所要求解的微分方程转化为相应的变分问题，即泛函求极值问题。
常见方法为把要分析的区域划分为很多三角形(每个三角形成为一个基元)，每个基元内的场用多项式来表达，然后加入不同基元间场的连续条件，就可得到整个横截面的场分布。
特点：较复杂---需要前处理(三角化，剖分)；后处理：（场分布，伪解剔除）（通用性强，精度高）根据该方法对于各种各样的电磁计算问题具有较强的适应能力性，所形成的代数方程矩阵求解容易、收敛性好。
主要缺点：对于形状和分布复杂的三维问题，由于其变量多和剖分要求细，往往因计算机内存而受到限制。
程序设计复杂、计算量较大。
适用范围：求解光波导的模式（有效折射率、色散、双折射、传输损耗等）。
现状：功能最强大的数值方法之一。
特别是上世纪90年代出现的矢量有限元方法，完全解决了有限元方法出现的伪解问题，大大降低了有限元法的后处理过程。
有限元光束传播法。
(3) 时域有限差分法时域有限差分法是近年来开始流行的一种数值模拟方法，它通过将麦克斯韦方程在时间空间上离散化的方法实现对电磁波传播的模拟。
它能够得到电磁波传输的瞬态（即时域）信息，通过傅里叶变换即可得到相应的频域信息。
时域有限差分法由K.S.Yee于1966年首先提出，此后经过众多学者的努力，使之不断完善，现已比较成熟。
但是，在许多方面它仍在继续发展，解决问题的能力和应用范围仍在不断地提高和扩大。
计算过程为：设置初始场，然后依时间步推进计算，并在每一时间步交替地计算每一离散点的电场和磁场。
特点：不需要矩阵运算，只需简单的加减乘除运算由前一时刻的场来获得下一时刻场的值。
而且，它还非常适合于并行计算，这正好与当今计算机的发展趋势相吻合，这就更加提高了时域有限差分法解决实际复杂问题的能力。
适用范围：计算光波导的模场分布、有效折射率；研究波导之间的连接、耦合问题。
注：主要用于一维和二维光波导的分析。
三维波导分析计算量稍大。
现状：ADI FDTD，可应用于各向异性介质，非线性介质，PML吸收边界(4)光束传播法(Beam Propagation Method，简写BPM)光束传播法是目前光波导器件研究与设计领域最流行的方法之一，其基本思想是在给定初始场的前提下，一步一步地计算出各个传播截面上的场。
光束传播法最早是由M.D.Feit等人于1978年研究光场及大气激光束传播时提出的。
最早的BPM是以快速傅里叶变换(Fast-Fourier Transform，称FFT)为数学手段实现的，称为FFT-BPM。
FFT-BPM源于标量波方程，只能得到标量场(即只能处理一个偏振分量)，不能分辨出场的不同偏振(TE模或TM模)以及场之间的耦合。
由于上述缺点，D.Yevick等人于1989年提出了一种新方法—有限差分光束传播法FD-BPM，用差分的方法将横截面上的场离散化。
这种方法已被成功地应用于分析Y型波导及S型弯曲波导中的光波传输，且对损耗的计算也得到了准确的结果；FD-BPM还被用于分析条形波导、三维弯曲波导、二阶非线性效应以及有源器件。
频域分析方面，同样可采用光束传播法进行分析：可采用相关函数法获得，还发展了一种称为虚轴光束传播法的方法，用于分析波导中的模式。
其实，BPM与FDTD有不少相似的地方。
其不同在于，FDTD每次都要同时计算整个波导的模场，而BPM只算一个面。
特点：计算量较小，应用范围非常广泛适用范围：计算光波导的模式、色散、双折射、传输损耗等；分析波导传输、连接、耦合，光栅的传输特性等。
4．数值方法发展趋势：方法融合现象明显(有限元法与光束传播法的结合形成了另外一种方法—有限元光束传播法（FE-BPM）。
)、相互推动(PML FDTD , BPM, FEM)。
第二节有限差分光束传播法基本原理光束传播法(BPM)的基本思想就是把波导沿着传播方向剖分成若干个截面，根据前一个或几个截面上的已知场分布得到下一个截面上的场分布。
BPM 理论来源于波动方程，波动方程是建立在Maxwell 方程基础上的。
Maxwell 方程的一般形式为(1a)()()0 t t t∂∇⨯+=∂B E(1b)()()()t t t t∂∇⨯=+∂D H J(1c)() t ρ∇⋅=D(1d)()0 t ∇⋅=B 式中，为电场强度，为磁场强度，为电矢量位移，为磁感应强度，E H D B 为电流密度矢量，为体电荷密度，为时间。
对于各向同性、非磁性、电J ρt 中性介质，有(2)0()(),()(),()()t t t t t t σμε===J E B H D E 式中，为电导率，为真空磁化率常数，为介电常数。
将式(2)代入式(1)σ0με有(3a)0()()0t t tμ∂∇⨯+=∂H E (3b)()()()t t t tσε∂∇⨯=+∂E H E 考虑到场对时间的依赖(4)()exp(),()exp()t i t t i t ωω==E E H H 式中，为复振幅，,x y z x y z E E E H H H =++=++Ei j k H i j k 为角频率，为单位虚数。
把式(4)代入式(3),有ωi(5a)00 i ωμ∇⨯+=E H(5b)()i σωε∇⨯=+H E 式(5a)可进一步写为( 6)0()()0i ωμ∇⨯∇⨯+∇⨯=E H 将式(5b)代入式(6),有(7)20()(0i σωμεω∇⨯∇⨯--=E E定义复相对介电常数(8)00(r i i i εσεεεεωε=+=- 将式(8)代入式(7)，就可以得到关于电场的矢量波方程(9)2()γ∇⨯∇⨯=E E 式中, 由下式表示γ (10)22220022()cωπγωμεεεελ=== 其中，为真空光速，为真空波长。
采用同样的过程，可以得到关于磁场矢c λ量波方程(11)21()()()εγε∇⨯∇⨯=∇⨯∇⨯+H H H 对于任何矢量，有G(12)2()()∇⨯∇⨯=∇∇⋅-∇G G G 又(13)0()0μ∇⋅=∇⋅=B H(14)() εεε∇⋅=∇⋅=⋅∇+∇⋅D E E E 从而可以进一步得到（15)0∇⋅=H (16)r r 11()()εεεε∇⋅=-∇⋅=-∇⋅E E E 把式(12)和式(16)代入到式(9)和(11)，可以得到(17)22r r 1()0εγε∇+∇∇⋅+=E E E(18)221()0εγε∇+∇⨯∇⨯+=H H H 考虑准TE 模()和准TM()模，有0z E ≈0z H ≈ (19a)221()0 r r εγε⊥⊥⊥⊥⊥∇+∇∇⋅+=E E E 2211()()0z zεεγεε⊥⊥⊥⊥⊥∂∂∇+∇⨯∇⨯+⨯⨯+=∂∂H H k k H H(19b)式中 (20)x y⊥∂∂∇=+∂∂i j(21) x y E E ⊥=+E i j(22)x y H H ⊥=+H i j 将式(19)写成分量形式如下2210r r x x y x r E E E E x x y εεγε⎡⎤⎛⎫∂∂∂∇+++=⎢⎥ ⎪∂∂∂⎝⎭⎣⎦(23a)2210 r r y x y y r E E E E y x y εεγε⎡⎤⎛⎫∂∂∂∇+++=⎢⎥ ⎪∂∂∂⎝⎭⎣⎦(23b) 22110y x x x x H H H H H y x y z z εεγεε∂⎛⎫∂∂∂∂∇+--+= ⎪∂∂∂∂∂⎝⎭ (23c)22110 y y x y y H H H H H x y xz zεεγεε∂∂⎛⎫∂∂∂∇+--+= ⎪∂∂∂∂∂⎝⎭(23d)进一步可以写为(24a)()()222222110 x x r x x r y r y r E E E E x x yz E E x y x yεγεεε⎡⎤∂∂∂∂+++⎢⎥∂∂∂∂⎣⎦∂⎡⎤∂∂+-=⎢⎥∂∂∂∂⎣⎦ (24b)()()222222110 y y r y y r x r x r E E E E y y xz E E y x y xεγεεε∂∂⎡⎤∂∂+++⎢⎥∂∂∂∂⎣⎦⎡⎤∂∂∂+-=⎢⎥∂∂∂∂⎣⎦ 222211x x x xH H H H y y x z z z εεεε⎡⎤⎛⎫∂∂∂∂∂∂++-+⎢⎥ ⎪∂∂∂∂∂∂⎝⎭⎣⎦ 2210yy x H H H y x y x γεε⎡⎤∂∂⎛⎫∂+-=⎢⎥ ⎪∂∂∂∂⎝⎭⎣⎦(24c)222211yy y y H H H H x x yz z z εεεε⎡⎤∂∂∂∂⎛⎫∂∂++-+⎢⎥ ⎪∂∂∂∂∂∂⎝⎭⎣⎦ 2210x x y H H H x y x y γεε⎡⎤⎛⎫∂∂∂+-=⎢⎥ ⎪∂∂∂∂⎝⎭⎣⎦ (24d)式(24)中，，，，是空间坐标，，的函数，x E y E x H y H x y ze 把，，，随的快速周期变化部分分离，令x E y E x H y H z (25a)(,,)exp() x x E x y z i z β=Φ- (25b)(,,)exp()y y E x y z i z β=Φ- (25c)(,,)exp() x x H x y z i z β=ψ- (25d)(,,)exp() y y H x y z i z β=ψ- (26)0 n cωβ=式中，为参考折射率，选择时应尽量接近导模的有效折射率，否则会影响计0n 算精度。