第一类边界问题的有限差分法探讨文档格式.docx
《第一类边界问题的有限差分法探讨文档格式.docx》由会员分享,可在线阅读,更多相关《第一类边界问题的有限差分法探讨文档格式.docx(16页珍藏版)》请在冰豆网上搜索。
三、有限差分法公式的推导:
把求解的区域划分成网格,把求解区域内连续的场分布用
网络节点上的离散的数值解来代替。
网格划分的充分细,才能
够达到足够的精度。
应用有限差分法计算静态场边值问题时,
需要把微分方程用差分方程替代。
用图形法解释如下:
y卒
有限差分网格划分
如图将区域D划分为正方形网格,网格线的交点称为节点,
两相邻平行网格线间的距离称为步距h。
然后,拉普拉斯方程
离散化,对于任一点0,有一阶偏导数:
X=Xo
(X。
h,y。
)-(X。
-h,y。
)
2h
而后,对于二阶偏导数:
x(x。
h/2,y。
)-x(x。
-h/2,y。
x=x。
h
(p-®
103
h2
对于丫轴同理:
因此拉普拉斯方程的差分格式为:
+<
P
紧邻边界节点的拉普拉斯方程的差分格式为:
亠亠
p(1p)q(1q)1p1q
其中p、q为小于1的正数;
1、2为边界上的节点,其值为对应边界
点处的值,是已知的。
具体如图:
4
应用数值计算解释(泰勒公式展开法):
1点电位的泰勒公式展开为
3点电位的泰勒公式展开为
1a
跨平■(訥方亦乔冰
当h很小时,忽略4阶以上的高次
项,得
同样可得
将上面两式相加得
口畀切丄护趴ip
F=~\——+—I=-——在上式中代入…'
'
一,得
肌=(F『+約+血+的+如/4
对于1,即F=0的区域,得到二维拉普拉斯方程的有限差分形式
旳=(幽+陷+伽十旳"
通过以上两种方法的推导,可得任意点的电位等于围绕它的四个点的电位的平均值。
当用网格将区域划分后,对每一个网络点写出类似的式子,就得到方程数与未知电位的网络点数相等的线性方程组。
已知
的边界条件在离散化后成为边界上节点的已知电位值。
四、差分方程组的解法方法一:
高斯一一赛德尔迭代法(简单迭代法)其步骤是先对每一网格点设一初值。
然后按一个固定顺序(一般点的顺序按“自然顺序”,即:
从左到右,从下到上)如图:
A
y
5
■
6
■■
7
■■■
8
9■
9
I_\
A—
之后,利用二维拉普拉斯方程的有限差分形式用围绕它的四个点的电位的平均值作为它的新值,当所有的点计算完后,用它们的新值代替
旧值,即完成了一次迭代。
然后再进行下一次迭代,如此循环。
如下
式:
i,r1,2,
:
(k1)_1:
•(k1)..:
(k1).「(k).「(k)
i,j4i-1,ji,j-1i1,ji,j1
(迭代公式1)
其式中的上角标(k)表示k次近似值,下脚标i,j表示节点所在位置,即第i行第j列的交点。
其中要特别注意:
在迭代过程中遇到边界点式,需将边界条件二f带入。
i,j"
j
循环迭代时当所有内节点满足以下条件时停止迭代:
(k1)一(k)
■■■■
i,ji,j
其中,W是预定的最大允许误差。
方法二:
逐次超松弛法
简单迭代法在解决问题时收敛速度比较慢,实用价值不大。
实际中常
采用逐次超松弛法(又称高斯一一赛德尔迭代法变形),相比之下它有两点重大的改进,第一是计算每一网格点时,把刚才计算得到的临近点的新值代入,即在计算(i,j)点的电位时,把它左边的点(i-1,j)和下面的点(i,j-1)的电位用刚才算过的新值代入,即:
(迭代公式2)第二,是引入“加速收敛因子”。
上式中的a即为“加速收敛因子”,且1《a<
2。
特别关注的是逐次超松弛法收敛的快慢与a有明显关系。
并且最佳a的取值随着条件的不同而不同,如何选择最佳a,是个复杂问题。
在计算时可以尝试求取最佳a值,以使计算快速。
五、应用计算机仿真有限差分法解决具体问题本次讨论我选择第四题为具体实例进行研究。
题目
如图所示,有一长方形的导体槽,a=20,b=15,设槽的长度为无限长,槽上有一块与槽绝缘的盖板,电位为100V,其他板电位为零,求槽内的电位分布。
b
①=100
►
ax
通过MATLAB进行仿真,运用有限差分法,源代码如下:
u=zeros(15,20);
i=2:
14;
u(i,20)=100;
forj=1:
20
fori=2:
14
u(i,j)=100/19*(j-1);
end
a=input('
pleaseinputa(1<
a<
2);
a二'
);
m=input('
pleaseinputm(1<
m);
m='
forn=1:
m;
forj=2:
19;
b=u(i,j);
c=u(i,j+1);
d=u(i+1,j);
e=u(i-1,j);
f=u(i,j-1);
g=0.25.*(c+d+e+f);
u(i,j)=b+a.*(g-b);
mesh(u);
通过对相关资料的查询及学习,我认为运用MATLAB编程就是要将有限差分法的基本公式实现。
以此为出发点可以让解析和编写过程相对简单化。
首先,先要对电场的边界条件进行设定:
使其满足题目的要求右侧边界电位为100,其余三边为0,由此也可以判断出此题目为第一类边界条件的问题。
同时,此处,也要特别注意,因为在MATLAB程序中,编写一个矩阵的默认顺序是从左到右,从上到下。
另外,命名矩阵时是先编写行,后编写列,因此,对于u(i,j)中,i代表的是第i行即对应纵坐标,j代表的是第j列即对应横坐标,与坐标表示有一定的差别,需要区分和注意。
之后是相对关键的编写部分,就是对电场内部的电位设定,此处,我通过分离变量法将u分解为u(x)和u(y),分别求解。
为让计算简便,我设定其内电场为线性变化的,即u(x)=kx+b,u(y)=0作为初值进行计算(应用分离变量法得到的结果是u(x)二BsinhmH/15u(y)=AsinmH/15。
可以通过此处在得出结果后进行对比和误差分析。
)又因为要在MATLAB^进行运算,所以首先要转化为i和j的形式,因此可得:
u(j)=kj+b和u(i)=0。
从而得到方程组:
u
(1)=k+b=0;
u(20)=20k+b=100。
解得:
u(j)=100/19*(j-1),即:
u(i,j)=100/19*(j-1)。
由此得到了内电场除去上下边界(衡为零)的电位方程,在MATLAB中我运用了
for循环语句来实现其内电场的赋值,源代码如下:
forj=1:
将电场各点的电位赋初值完成后,就需要带入公式进行计算了。
在之
确定。
关于加速因子的最佳值确定问题,我通过学习和查阅相关资料
得到了两种方法,如下:
一是逐步搜索法;
将a的取值区间(1,2)进行M等分,a分别取1+1/M,1+2/M,,
1+(M-1)/M,通过上面提到的迭代公式2依次对同一个精度要求求出迭代次数k的值,并从中选出“加速收敛因子”a的最佳值,具体
步骤如下:
步骤1给定等分数M和精度要求£
的值,令a的初始值为1;
步骤2令p=1,2,3,,,M-1,重复步骤3-5;
步骤3ap=1+p/M;
步骤4按照如下公式迭代
ot
.:
(k1)=.:
(k)._■:
(k1).「(k1).「(k).「(k)_4(k)
i,ji,ji-1,ji,j-1i1,ji,j14i,j
找出符合精度要求£
的迭代次数kp;
步骤5比较找出kp值最小的ap为最佳的加速收敛因子a的值。
二是黄金分割法:
依据黄金分割法的思想,通过计算机自动选取最佳加速收敛因子
a的近似值,具体步骤如下:
步骤1对(1,2)区间第一次0.618的分割,区间分界a1=1,b1=2,在(a1,b)区间分割出黄金点point1=a1+0.618(b1-a1),进行迭代公式2的迭代,求出迭代次数K值,如果迭代次数没有超出规定的发散
常数,迭代结束,否则转为步骤2
步骤2在(1,1.618)和(1.618,2)之间进行第二次的黄金分割,找出分割点point2=a2+0.618(b2-a2),其中a?
和b?
是新分割区间的左右边界。
找出迭代次数最少的a。
以上两种方法是比较具有实践和研究价值进行的,在这里介绍和分享给大家探讨。
在我分析过后认为第二种较好,因为其可以通过计算机自动选取到最佳的加速收敛因子,但是实施起来比较复杂,要通过MATLAB或其它语言进行编码实验得到,我初步试验过,没能成功,就运用了第一种方法进行了计算,但是,我发现其运算过于复杂,笔算正确的成功率比较低,而运用电脑编程要重复多次运算,虽然,计算机运算快速且方便,但是,这样的源代码同样过于复杂,可以进一步进行深入的研究。
本次的题目解答中,由于是应用计算机计算,所以计算速度和计算效率在不同的a取值中没有明显区别。
因此,我就只是设定了一个1~2之间的一个数字:
1.8为“加速收敛因子”a的值了,进行了近似计算,实验误差在可接受的范围。
然后,是对于迭代次数的确定。
当迭代满足停止条件:
w(k*1)—申(k)£
\w
i,ji,jVV
时,对应的迭代次数即为最终的迭代次数,这里的w为设定好的允许
误差,为从简方式,我选择运用输入迭代次数,并通过多次输入对比输出结果选取最佳迭代次数的方法来选取,存在了一定的误差。
但是,相比于我运用上面的公式求解计算迭代次数,进而得出的图形的误差比要小许多,我一直尝试进行改进,但一直没能很好的解决误差过大的问题:
最主要的问题是边界值一直不满足我的设定值。
希望通过进一步的研讨和大家的交流能过得到最优方式。
通过以上的分析,我得出了进一步的编程代码:
a=input('
pleaseinputa(1<
a='
m=input('
pleaseinputm(1<
其中a和m分别为最佳“加速收敛因子”和迭代次数,均为我们输入数值计算的方式。
先要预算出一些取值范围,然后通过分别输入,对比输出结果进而来选取最优解。
之后的源代码的编写思想就是通过迭代计算公式来实现对电场内部各点电位值的计算,因此,在输入“加速收敛因子”和“迭代次数”的数值后,我通过for的循环语句实现了对电场内部的电位迭代求解,并用图形形式表示出来。
其源代码如下:
forn=1:
u(i,j)二b+a.*(g-b);
endend
上述代码中,b,c,d,e,f分别场内一点及其左右上下各点的电位值,之后的g和u(i,j)则实现了逐次超松弛法的迭代,最后用mesh实现图形的输出。
六、结果和误差的对比分析
由分离变量法得到的初步结果是u(x)二Bsinhmn/15和u(y)=AsinmQ/15,从而应有u=Csinh(mn/15)sinm(n/15)的结果,从其表达式的形式,我们可以简单推断出电场的电位在x轴上的分布是成虚指数函数的形式呈现,而在y轴上是成正弦函数呈现的,这和通过有限差分法的计算机仿真的图形基本符合,如下图:
因此,可以基本确定通过计算机仿真得到的结果的正确性。
另外,可能存在较大误差的方面是跌代次数,因为迭代次数是我们人工输入的,所以判断其正确与否要通过输出的结果来判定、来选择,下面几幅图是当“加速收敛因子”a=1.8时,不同的迭代次数m的对应图形:
m=1时:
Qif
2D
m=10时:
m=20
m=30
m=32
ifl.7
m=40
通过上面迭代次数m取不同的值对应的不同图形分析可以看出,在m=1,10,20时期计算结果都存在较大误差,图形与分离变量法得到的结果相差很大,有明显的不规则部分,不成现光滑的曲线。
而当m=30后图形与分离变量法计算到的图形近似,而且图形规则,曲线光滑,而且随着m值的增大图形的变化非常微小了,因此,我选择m=32
为迭代次数,所得到的结果误差比较小,结果相对准确。
七、研究结论
学习了运用有限差分法计算第一类边界条件的电位值,通过计算
机仿真进行了实际例题的求解,通过与分离变量法求解的对比,证明可行性和有效性,适用于解决第一类边界问题,而且更加方便和快捷。
同时,“加速收敛因子”最佳解问题比较复杂,不过通过计算机仿真的过程中,不同数值“加速收敛因子”的加速效果不易体现出来。
跌代次数对计算结果又较大影响,迭代次数越多,结果越准确,误差越小,但相应的计算也会更加复杂,因此只要迭代次数得到得结果满足一定的允许误差范围即可。
附录:
1三类边界条件:
第一类:
已知电位在边界上的数值。
第二类:
已知电位的导数在边界上的数值。
第三类:
已知电位在部分边界上的数值和在其
余边界上的导数的数值。
2解析法包括:
分离变量法,镜像法,复变函数法,保角变换法,格林函数法。
3数值法包括:
数值积分法,有限差分法,有限元法,矩量法。
4解析法的优点:
由有限个项构成闭合解或无穷级数,通常具有鲜明的物理意义。
缺点:
解题范围窄小,对边界形状十分挑剔。
数值法的优点:
解题范围宽广,对边界的形状没有限制。
缺点:
数据离散,物理意义深藏其中。
心得体会:
有限差分法的公式套用很简单,记住结果很简单,但是要想深刻理解其内容和灵活应用则比较困难,需要一定时间的深刻学习和做题应用,同时,这个方法确实是一个非常方便和有效地解决三类边界条件问题的工具。
MATLAB在电磁场方面的应用相当广泛,其内拥有的很多函数都可以模拟静电场中的实例,就比如这次的研讨就应用了矩阵,通过有限差分法来模拟静电场的电位分布,相信将来还会有更多的有意思,有研究价值计算和模拟电磁场的相关知识。
是学好电磁场必须要牢固掌握的工具。