电磁场数值计算7-西安交通大学电气工程学院

1、 第七章 有限差分法 差分法是运用最早的一种数值计算方法。因为它的简便、直观、易编程等优点而一直被广泛应用，它可以求解均匀或不均匀、线性或非线性媒质中的位场，可以求解恒定场或似稳场。也可以用于瞬态场中的时间差分。§71 差分运算的基本概念 7.1.1 一阶差分与差商 设函数，自变量x有一微小增量，则的相应增量是 (7-1) (7-2)式(7-1)称为的一阶差分，它与微分不同，由于是有限量的差，故又称为有限差分。式(7-2)称为的一阶差商，或有限差商。 向前差分 向后差分 中心差分相应的有 一阶向前差商 (7-3) 一阶向后差商 (7-4) 一阶中心差商 (7-5)只要h足够小，那么，

2、差商和导数、差分与微分之间的差异将很小，它们之间的关系为 一阶微分 (7-6) 一阶导数 (7-7) 哪一种差商误差最小呢？通过泰勒公式，函数及在x点展开级数，便可以得知： (7-8) (7-9)可见，式(7-5)和(7-6)都截断于项，略去项以及h的更高幂次的项。而将式(7-8)式(7-9)，得到一阶中心差分式(7-7)是 （7-10）截断于项，略去项以及h的更高幂次的项。因此，中心差商的截断误差最小，具有二阶的精确度。7.1.2 二阶差分与差商二阶导数可以从式(7-8)+(7-9)得到，即截断于项，略去项以及h的更高幂次的项，得到 （7-11）也可以写成同理可以得到 （7-12） （7-1

3、3）可以看出，这是中心差商，具有较高的精度。 将二阶导数的近似表达式代入泊松方程、波动方程等，就将原微分方程转化为代数方程。7.1.3 离散单元 离散单元的形状多为图71 差分法离散单元形状§72 泊松方程的差分形式7.2.1 建立差分格式1、 平行平面场 泊松方程 (7-14)图7-2 五点差分格式如图7-2所示，按中心差分格式可以导出（设）代入泊松方程，得到用节点坐标写出“五点差分格式” (7-15)若 (7-16)2、 轴对称场泊松方程 (7-17)设，令轴线处开始，任一点的半径为， 当，差分格式为（作业7-1：推导轴对称场差分格式） (7-18)当，点落在（z轴）上，这时方程

4、（7-17）中等于零，出现第二项为无限大的现 象，为保证方程每一项有限，必须满足z12340对称线由洛必塔法则得知因此，当时，对称轴上的微分方程为 相应的差分格式为 图7-3 轴对称差分格式 (7-19)7.2.2 媒质分界面边界条件的离散格式分界面上面电荷为零，因此有 应用高斯散度定理是S面的外法线方向，即垂直于每一周边l向外的单位矢量，展开后 41230 图7-4 媒质分界面差分格式代入高斯散度定理，得到 (7-20)结点2、4位于分界面上，可以理解为在分界面上取等效介电常数。 上式可以写为（） 同理，可以分析介质角点处图7-5的差分格式 (7-21)以及分界面对于网格呈对角线形态，这时需

5、要引入两个辅助点，如图7-6所示，可以推得 (7-22)式中，。 M12340N 图7-5 介质角点处的差分格式图 图7-6 分界面对于网格呈对角线形态7.2.3 第一类边界条件的离散格式 图 7-7 第一类边界条件 图7-6 分界面对于网格呈对角线形态对于第一类边界条件， 分两种情况： （1）若网格结点落在边界上，直接赋值即可。（2）若网格结点不落在边界上，有两种处理方法：一种是直接转移法，取最靠近O点的边界结点上的函数值作为O点的函数值，即如图7-7所示另一种采用线性插值法，先判断x、y方向上哪一个边界点最靠近O点。设结点1最靠近结点O，则可用插值公式得到 (7-23) 若结点2最靠近结点

6、O，则可用插值公式得到 (7-24) 最精确的算法：采用不同步长的五点差分格式 (7-25)其中，7.2.4 第二、三类边界条件的离散格式（1）若边界条件满足 通常表示对称边界。当边界线与网格线重合时，根据对称性，边界外与结点1对称位置的虚拟结点4的电位满足，即强迫。因此，当步长相等时，有 (7-26)(2) 第三类边界条件 图7-8 第三类边界条件当时为第二类边界条件。对于比较复杂的边界，边界的法线方向既不平行于、也不垂直于坐标轴，如图7-8所示，可以用以下简便方法： 过O点作边界的垂线QP，（P、Q均为交点），对O点有 (7-27)P点不是结点，用结点1、2上的插值给出，由得到 （7-28

7、）代入式(7-27)，得到其中，。又因为 所以 代入边界条件得到结点o的差分格式为 (7-29)§73 差分方程的求解 方程的特点：（1）系数矩阵是稀疏矩阵，每行非零元素不超过5；（2）稀疏矩阵是对称、正定矩阵（但是，当边界与网格结点不重和时，对称性被破坏）可以采用直接法和迭代法。通常方程阶数小于500，可以用直接法，否则用迭代法。迭代法有：高斯迭代、高斯赛德尔迭代、超松弛迭代等。1、高斯迭代 (7-30)2、高斯赛德尔迭代 (7-31)3、超松弛迭代 (7-32)式中，。若，就是高斯赛德尔迭代。 例题 用有限差分法求解矩形金属波导中的电磁波分布和截止波长。 这是稳态时变场问题。设：

8、波导壁由理想导体构成；且波导内介质均匀、线性、各向同性；远离激励源；无反射。 波导中传播的TE波（传播方向上无电场分量）或TM波（传播方向上无磁场分量），无论哪种情况，可归结为求相应的（TE波，求得后，根据麦克思韦方程可求得）或（TM波，求得后，根据麦克思韦方程可求得）。 用表示场的纵向分量或，则波导内的波动方程 (7-33)式中，未知（？）。由于波导是理想导体，沿面电场的切向分量为零，只有法向分量，即：在与x轴平行的波导壁上，有，在与y轴平行的波导壁上，有。由于因此，TE波： (7-34) TM波： (7-35)结点在波导内部的迭代格式： (7-36)结点在边界上，若是TE波，根据前面所述，

9、有 (7-37)若是TM波，则有 (7-38)称为特征值，未知。是截止波长。这样，时变场中的波动方程的特征值问题近似的变为离散系统中的代数特征值问题。 由于和都不知道，因此要用双重迭代的方法求解，先给和赋初值，用超松弛方法 (7-39)迭代收敛后，再对进行计算 (7-40)代表以第个结点为中心所占的面积。如图所示， 对于内结点 不是角点的边界点 对角点的边界点 式(7-38)的推导过程：给波动方程乘以，再在全域上作面积分解得 从式(7-35)中也能求，这是由局部点求得，全部点求完后再取其平均值，但这种方法误差较大。而(7-39)是对全域进行面积分后取平均值，计算精度相对较高。 求得后，继续应用

10、式(7-38)，直到和收敛到应有的精度。（第十三次课）§74 标量时域有限差分法（参考吕英华计算电磁学的数值方法） 有限差分法既能用于静态和时谐电磁场问题，也能用于时域电磁场问题。本节介绍抛物型和双曲型偏微分方程的有限差分法。7.4.1 抛物型方程的有限差分法 抛物型方程是描述传导、扩散等非平衡过程的，定解问题是 （7-41） 所有x位置上的初始值 （7-42） 所有t位置上的边界值 （7-43）式中都是已知的连续函数，并满足相容性条件：即当时的边界值等于边界的初始值。 如图7-8所示，设h为空间步长，为时间步长，那么 j=0,1,2,J ， i=0,1,2,I一、显示差分格式仍然用

11、泰勒级数在点的展开公式，略去截断误差项， 一阶向前差商 二阶中心差商 代入式（7-40），得到 整理后 （7-44）式中。可以看出，它与泊松方程的五点差分格式不同，第层上的任意节点的值可以由第i层的三个相邻节点处的值来决定。解的误差为。 初始条件离散化（除了两个边界点以外的所有场点） （7-45）边界条件离散化（包含t=0时的边界点） （7-46）二、隐式差分格式 从算法的稳定性考虑，显示差分格式要求h和满足一定的关系，通常较小，因此计算量大。隐式差分格式可以有效地解决这一问题。 对时间采取向前差商，泰勒级数在点，即点展开，略去截断误差代入式（7-40），得到 （7-47）整理得（） （7-4

12、8）解的误差为，这是关于层上未知量的联立方程组，不再是显示格式，需要求解方程组，故称为隐式格式，其最大优点是无条件稳定。三、六点差分格式 将显示格式与隐式格式结合起来，可以构成精度更高的格式。设参量，用式（7-47）+式（7-44），得到整理得到 （7-49）方程表达了i和i+1两层上相邻6个节点的函数关系，故称为6点差分格式，可以证明，它的算法误差为，比前两种算法的时间精度提高了一个数量级。可以看出，当时，是显示差分格式，时，是隐式差分格式，需要联立求解方程组。通常取7.4.2 双曲型方程的有限差分法 双曲型方程描述电磁波的传播，定解问题一般为 （7-50） 所有x位置上的初始值 （7-51

13、） 所有t位置上的边界值 （7-52）式中都是已知的连续函数，并满足相容性条件：取矩形网格，在处，二阶中心差商为代入式（7-50），略去截断误差，得到 一、显示格式 设，同上推导，差分格式为 （7-53） （7-54） （7-55）（第一种：由于 ，用一阶差分近似，在t = 0的差分为 上式可以作为一阶导数初始值。第二种：也可以从式（7-53）得到，当t =0，即i =0时 第二种初速度条件离散时采用了二阶差商格式，所以比第一种精确。） (7-56)二、隐式格式把时间第i-1层、i层、i+1层的中心差商加权平均将上两式代入式（7-50），略去截断误差，整理可得当时，是显示差分格式，时，是隐式差分格式，需要联立求解方程组。通常取。可以证明，隐式格式是无条件稳定的，精度。初始条件与边界条件不变 以上都是在一维场中的差分格式，若是二维、三维场问题，差分方法相同，只是多增加了y、z分量，只需要将含有x的所有项分别用y、z替代，叠加后即可得到表达式更加冗长的公式。7.4.3 稳定性分析 上述方法的稳定性与时间步长、空间步长有关。一般的稳定性条件为 三维 若步长不相同 二维 7.4.5 网格色散误差 既使算