武汉大学工程电磁场及高电压综合实验一、题目有一极长的方形金属槽，边宽为1cm，除顶盖电位为100sinπxV外，其他三面的电位均为零，试用差分法求槽内电位的分布。

二、解题原理：均匀媒质中的有限差分法我们在求解场的分布时，当边界形状比较复杂时，解析分析法不再适合了，我们可以采用数值计算的方法，数值计算法的基本思想，是将整体连续的场域划分为若干个细小区域，一般称之为网格或单元，如图1所示，然后用所求的网格交点（一般称为节点或离散点）的数值解，来代替整个场域的真实解。

因而数值解，即是所求场域离散点的解。

虽然数值解是一种近似解法，但当划分的网格或单元愈密时，离散点数目也愈多，近似解（数值解）也就愈逼近于真实值。

实解。

在此处键入公式。

图1场域的剖分，网格节点及步长(一）、场域的剖分、网格节点及步长由边界Γ所界定的二维平行平面场（见图1），若采用直角坐标系则可令该场处在xoy 平面内。

所谓场域的剖分就是场域的离散化，即将场域剖分为若干个网格或单元。

最常见最简单的剖分为正方形剖分，这种剖分就是在xy 平面上作许多分别与x 轴及y 轴平行的直线，称为网格线。

网格线的交点称为节点或离散点，场域内的节点称为内节点，场域边界上的节点称为边界节点。

两相邻网格线间距离称为步长，一般以h 表示。

若步长相等则整个场域就被剖分为许多正方形网格，这就是正方形剖分。

节点（离散点）的布局不一定采用正方形剖分，矩形剖分也常采用，正三角形剖分偶尔也被应用，不过最常见的最简单的仍然是正方形剖分。

（二）、差分与微分从前面的分析可知，稳恒电、磁场的求解问题，归根到底是求解满足给定边界条件的偏微分方程（泊松方程或拉普拉斯方程）的解的问题所谓差分方法，就是用差商近似代替偏微商，或者说用差分代替微分，从而把偏微分方程转换为差分方程，后者实际上为代数方程。

因此这种转化有利于方程的求解。

下面分别对一阶及二阶的差分公式进行推导。

首先回顾有关偏导数的定义，有00(,)(,)(,)(,)lim lim x x f f x x y f x y f x y f x x y x x x→→∂+---==∂ (1) 因此当|x| 充分小时，可近似地用(,)(,)f x x y f x y x +- 或(,)(,)f x y f x x y x-- 代替fx∂∂，所谓差分公式，即是基于上述观点推得的。

设图1所示场域中的位函数为A ，任取一网格节点0，它在xy 平面上的坐标为（x ,i i y ），记节点0的矢量磁位为,i j A ，并把与节点0相邻的其他四个节点1、2、3、4的矢量磁位分别记为1,i j A +、,1i j A +、1,i j A -、,1i j A -，将节点0处函数A 的一阶偏微商Ax ∂∂，用1、0两点函数值的差商1,,i j i j A A h+-近似代替，则有1,,0()i j i j A A Ax h+-∂≈∂（2） 式（2）中之差商，称为向前差商。

上述一阶偏微商也可用0、3两点函数值的差商,1,i j i jA A h--近似代替，称为向后差商，得,1,0()i j i j A A Ax h--∂≈∂（3） 同理，对于偏微商0()Ay∂∂也可分别用向前或向后差商近似代替，所得结果为 ,1,0()i j i j A A Ay h+-∂≈∂(4) 或 ,,10()i ji j A A Ay h--∂≈∂(5) 式（2）~式（5）就是一阶差分公式。

二阶差分公式可以在一阶差分的基础上进一步推出，0点处的二阶偏微商2002()()A Ax x x∂∂∂=∂∂∂，如对一阶差商再取差商则得二阶差商为 1,,,1,i j i ji j i jA A A A hhh+----用上式近似代替二阶偏微商就是要推演的二阶差分公式。

又由于1,,,1,i j i ji j i jA A A A hhh+----=1,,1,22i j i j i jA A A h +--+（6）所以21,,1,0222()i j i j i j A A A A x h +--+∂≈∂（7）同理2,1,,10222()i j i j i j A A A A y h +--+∂≈∂（8）式（7）及式（8）即为二阶差分公式。

需要说明的是，用差分近似代替偏微分，必定会产生误差。

理论分析表明：其误差与步长2h 成比例，因此若网格剖分得愈小，步长h 就愈小，从而引起的误差也愈小。

(三)、均匀媒质中泊松与拉普拉斯方程的差分离散格式设图2所示的平行平面场，场域每边长均为b,场域内电流密度为δ，媒质磁导率为0μ，边界上的矢量磁位值已知，求域内矢量磁位。

所提问题为第一类边值问题，则20A μδ∇=-|()(5,6,16),i A f s i Γ==（9）应用差分法解图2所示场域的步骤如下：图2 边界与网格线重合第一步，若采用正方形网格剖分，即将场域剖分为具有四个内点（即点1,2,3,4），边界与网格线重合的九个离散单元（九个网格）。

这就是把连续的场域进行离散化，从而将求解场域内矢量磁位函数的问题，转化为求1,2,3,4各内点的矢量磁位值的问题。

第二步，根据式（7）、式（8）的二阶差分公式，列出1,2,3,4各内点泊松方程的差分表达式。

因为h=b/3，则各内点1,2,3,4的矢量磁位1234,,,A A A A 的差分方程为内点1： 223166104A A A A A h μδ+++-=-（10） 内点2： 21497204A A A A A h μδ+++-=-（11） 内点3： 2115134304A A A A A h μδ+++-=-（12） 内点4： 2312102404A A A A A h μδ+++-=-（13） 由于边界各点的矢量磁位值为已知，即16A =16f ，6699771515131312121010,,,,,,A f A f A f A f A f A f A f =======。

将上述已知量代入式（10）~式（13）之中，并将各已知量移至等式右端，则得泊松方程的差分离散格式。

即14A - +2A +3A =20616h f f μδ---1A 24A -4A + =2079h f f μδ---（14） 1A 34A -4A + =201315h f f μδ--- 2A +3A 44A -=201012h f f μδ---拉普拉斯方程的差分离散格式则更趋简单，从式（14）所列的矢量磁位的代数方程可以看出，待求量（1234,,,A A A A ）的个数与方程式的数目一致。

各方程式左边为待求量，右边各项是乘积20h μδ与边界节点上矢量磁位之值的代数和，它们均是已知的。

由于待求量的数目与方程数目一致，且各方程右端项不尽为零，故此代数方程组有非零解。

第三步，解代数方程组。

当内点较少时，可直接用待元消去法或列式法、张弛法等进行手算；当内点较多时，即内点数不是几个，十几个，而是成百个，上千个时，手算几乎不可能，这就必须借助计算机进行计算，求解高阶方程组的方法有赛德尔迭代法及超松弛代法等等。

我们运用分离变量法求得其解析解，若用差分法则可直接求得场域中离散点上电位的近似值。

首先对场域进行等距剖分，例如，步长h=0.25，对于正方形场域则可使用网格线自边界处起始，平行于y 轴的网格线x=i h （i=0,4;j=0,4）由边界条件给出，其内部节点的电位值,i j ϕ（i=1,2,3;j=1,2,3）则待求。

电位函数所满足的拉普拉斯方程的差分离散格式为1,,11,,1,4i j i j i j i j i j ϕϕϕϕϕ++--+++=即 ,1,,11,,11()4i j i j i j i j i j ϕϕϕϕϕ++--=+++（15）对于本题的网格剖分，i,j=1,2,3，则式（5-54）即为待求的内部节点上的电位值所应满足的代数方程组。

题设的边界条件：0,4,,0,4,100sin 4j j i i i πϕϕϕϕ===，代入相应的代数方程之中。

如i,j=1，上式即为1,12,11,20,11,01()4ϕϕϕϕϕ=+++其中边界值0,1ϕ=1,0ϕ=0应代入方程之中，而2,1ϕ与1,2ϕ则为待求量。

求解代数方程组得,i j ϕ（i=1,2,3;j=1,2,3），此即电场中电位分布的数值解。

解代数方程组方法较多，若采用赛德尔迭代法，则可将式（15）改写为(1)()()(1)(1),1,,11,,11()4n n n n n i j i j i j i j i j ϕϕϕϕϕ+++++--=+++（16）的形式，式中标号（n ）为第n 次计算值，（n+1）为第n+1次的计算值。

运用式（16）时，可从j=1开始，依次对i=1,2,3进行计算；再对j=2,i=1,2,3进行计算；最后当j=3时，对i=1,2,3进行计算。

每完成一次对i 或j 的循环，(),n i j ϕ全部换为(1),n i j ϕ+，这叫做完成一次迭代。

经过十数次或数十次这样的迭代，当两次邻近的迭代值相差足够小时，则可认为得到了电位函数的近似数值解。

由于计算格式十分有规则，因此上述步骤实际上往往在计算机上进行，这时取步长h 为更小值，可提高数值解的精度。

三、编写MATLAB 程序按照以上分析的结果，我们编写MATLAB 程序来计算，并对网格划分程度不同的情况进行比较。

1 划分场域的网格为5行5列：首先，我们取步长较小，构造一个五行五列网格来计算场域的分布： 程序：改变不同的迭代次数值t，我们得到了不同迭代次数下的计算结果：线性赋初值的结果：迭代1次的结果：迭代5次的结果：迭代10次的结果：迭代20次的结果：迭代50次的结果：迭代100次的结果:通过比较，我们发现，在迭代次数到达20次时，其误差已经很小，继续迭代计算结果，发现迭代50次和迭代100次后的结果已经没有什么差别（只考虑到小数点后四位）。

2 划分场域的网格为11行11列增加步长，我们把原区域划分为一个11行，11列的网格，计算程序如下：改变不同的迭代次数值t，我们得到了不同迭代次数下的计算结果：线性赋初值的结果：1次迭代后的计算结果：5次迭代后的计算结果；10次迭代后的计算结果：50次迭代的结果：70次迭代的结果：100次迭代的结果：150次迭代后的结果：200次迭代的结果：我们对比结果可以知道，在迭代到150次和迭代到200次后，其结果已经一样（同样只考虑到小数点后四位）。

两次结果对比：对比不同网格下的电位分布情况，我们能够清楚的看见，当划分的网格越小时，既步长越短时，我们得到的电位分布情况就越接近实际的分布情况。

这样，如果我们继续缩小我们的步长，进一步地把我们的网格划小，并通过计算来得到我们就可以得到场域分布情况。