有限差分法实验报告

工程电磁场实验报告——有限差分法用超松弛迭代法求解 接地金属槽内电位的分布一、实验要求按对称场差分格式求解电位的分布 已知：给定边值：如图1-7示 图1-7接地金属槽内半场域的网格=ϕ= V100 ϕ 0=ϕ0=ϕ给定初值)()(.1j 401001j p12j i -=--=ϕϕϕ 误范围差: 510-=ε计算：迭代次数N ，j i ,ϕ，将计算结果保存到文件中二、实验思想有限差分法有限差分法（Finite Differential Method ）是基于差分原理的一种数值计算法。

其基本思想：将场域离散为许多小网格，应用差分原理，将求解连续函数ϕ的泊松方程的问题转换为求解网格节点上ϕ的差分方程组的问题。

泊松方程的五点差分格式)(414243210204321Fh Fh -+++=⇒=-+++ϕϕϕϕϕϕϕϕϕϕ当场域中,0=ρ得到拉普拉斯方程的五点差分格式)(41044321004321ϕϕϕϕϕϕϕϕϕϕ+++=⇒=-+++差分方程组的求解方法(1) 高斯——赛德尔迭代法][)(,)(,)(,)(,)(,2k 1j i k j 1i 1k 1j i 1k j 1i 1k j i Fh 41-+++=+++-+-+ϕϕϕϕϕ (1-14)式中：⋅⋅⋅⋅⋅⋅=⋅⋅⋅⋅⋅⋅=,2,1,0,2,1,k j i ，迭代顺序可按先行后列，或先列后行进行。

迭代过程遇到边界节点时，代入边界值或边界差分 格式，直到所有节点电位满足εϕϕ<-+)(,)(,k j i l k j i 为止。

图1-4 高斯——赛德尔迭代法（2）超松弛迭代法][)(,)(,)(,)(,)(,)(,)(,k j i 2k 1j i k j 1i 1k 1j i 1k j 1i k j i 1k j i 4Fh 4ϕϕϕϕϕαϕϕ--++++=+++-+-+ (1-15)式中：α——加速收敛因子)21(<<α 可见：迭代收敛的速度与α有明显关系三、程序源代码#include<iostream.h> #include<math.h> #include<iomanip.h> double A[5][5]; void main(void) { double BJ[5][5];//数组B 用于比较电势 int s[100];//用于储存迭代次数double d[100];//用于记录所有的加速因子 d[0]=1.0;int i,j,N=0,M=0,x; for(i=0;i<100;i++)d[i]=0.01*i+d[0];//加速因子从1.0到2.0之间的20个数！double w[100][10];int P,Q;for(P=0;P<4;P++)for(Q=0;Q<5;Q++)A[P][Q]=0;for(P=0;P<5;P++)A[4][P]=100;cout<<"数组A的所有元素是:"<<endl; for(i=0;i<5;i++)for(j=0;j<5;j++){cout<<A[i][j]<<setw(6);if((5*i+j+1)%5==0)cout<<'\n';}int pp=0;for(x=0;x<100;x++){do{for(i=0;i<5;i++)for(j=0;j<5;j++)BJ[i][j]=A[i][j];for(i=1;i<4;i++)for(j=1;j<4;j++)A[i][j]=BJ[i][j]+(d[x]/4)*(BJ[i+1][j]+BJ[i][j+1]+A[i-1][j]+A[i][j-1]-4*BJ[i][ j]);//迭代公式for(i=1;i<4;i++){for(j=1;j<4;j++)if(fabs(A[i][j]-BJ[i][j])<1e-5)pp++;}N++;}while(pp<=9);pp=0;for(i=0;i<3;i++)w[M][i+1]=A[1][i+1];for(i=3;i<6;i++)w[M][i+1]=A[2][i-2];for(i=6;i<9;i++)w[M][i+1]=A[3][i-5];s[M]=N;M++;N=0;int P,Q;for(P=0;P<4;P++)for(Q=0;Q<5;Q++)A[P][Q]=0;for(P=0;P<5;P++)A[4][P]=100;}int min=s[0];int p,q;cout<<"输出所有的加速因子的迭代次数:"<<'\n';for(q=1;q<100;q++){// cout<<s[q]<<setw(6);// if(q%12==0)// cout<<'\n';if(min>s[q]){min=s[q];p=q;}}cout<<endl;if(min==s[0])p=0;cout<<"最佳加速因子a=";cout<<d[p]<<'\n';cout<<"迭代次数为："<<min<<'\n';cout<<"最佳收敛因子对应的各个格内点的电位为："<<'\n';for( i=1;i<10;i++){cout<<w[p][i]<<'\t';if(i%3==0)cout<<'\n';}cout<<'\n';}四、程序框图迭代解程序框图五、结果分析迭代收敛的速度与α的关系收敛因子（α）1.0 1.7 1.8 1.83 1.85 1.87 1.902.0 迭代次数（N ）>1000269174143122133171发散最佳收敛因子的经验公式：)sin(p120πα+=（正方形场域、正方形网格）220q 1p 122+-=πα（矩形场域、正方形网格） 程序执行结果如下。

工程电磁场实验报告——有限差分法用超松弛迭代法求解 接地金属槽内电位的分布一、实验要求按对称场差分格式求解电位的分布 已知：给定边值：如图1-7示 图1-7接地金属槽内半场域的网格 给定初值)()(.1j 401001j p12j i -=--=ϕϕϕ 误范围差: 510-=ε计算：迭代次数N ，j i ,ϕ，将计算结果保存到文件中二、实验思想有限差分法有限差分法（Finite Differential Method ）是基于差分原理的一种数值计算法。

其基本思想：将场域离散为许多小网格，应用差分原理，将求解连续函数ϕ的泊松方程的问题转换为求解网格节点上ϕ=ϕ V100 0=ϕ0=ϕ的差分方程组的问题。

泊松方程的五点差分格式)(414243210204321Fh Fh -+++=⇒=-+++ϕϕϕϕϕϕϕϕϕϕ当场域中,0=ρ得到拉普拉斯方程的五点差分格式)(41044321004321ϕϕϕϕϕϕϕϕϕϕ+++=⇒=-+++差分方程组的求解方法(1) 高斯——赛德尔迭代法][)(,)(,)(,)(,)(,2k 1j i k j 1i 1k 1j i 1k j 1i 1k j i Fh 41-+++=+++-+-+ϕϕϕϕϕ (1-14)式中：⋅⋅⋅⋅⋅⋅=⋅⋅⋅⋅⋅⋅=,2,1,0,2,1,k j i ，迭代顺序可按先行后列，或先列后行进行。

迭代过程遇到边界节点时，代入边界值或边界差分 格式，直到所有节点电位满足εϕϕ<-+)(,)(,k j i l k j i 为止。

（2）超松弛迭代法][)(,)(,)(,)(,)(,)(,)(,k j i 2k 1j i k j 1i 1k 1j i 1k j 1i k j i 1k j i 4Fh 4ϕϕϕϕϕαϕϕ--++++=+++-+-+ (1-15)式中：α——加速收敛因子)21(<<α 可见：迭代收敛的速度与α有明显关系三、程序源代码#include<> #include<> #include<> double A[5][5]; void main(void)图1-4 高斯——赛德尔迭代法{double BJ[5][5];//数组B用于比较电势int s[100];//用于储存迭代次数double d[100];//用于记录所有的加速因子d[0]=;int i,j,N=0,M=0,x;for(i=0;i<100;i++)d[i]=*i+d[0];//加速因子从到之间的20个数！double w[100][10];int P,Q;for(P=0;P<4;P++)for(Q=0;Q<5;Q++)A[P][Q]=0;for(P=0;P<5;P++)A[4][P]=100;cout<<"数组A的所有元素是:"<<endl;for(i=0;i<5;i++)for(j=0;j<5;j++){cout<<A[i][j]<<setw(6);if((5*i+j+1)%5==0)cout<<'\n';}int pp=0;for(x=0;x<100;x++){{for(i=0;i<5;i++)for(j=0;j<5;j++)BJ[i][j]=A[i][j];for(i=1;i<4;i++)for(j=1;j<4;j++)A[i][j]=BJ[i][j]+(d[x]/4)*(BJ[i+1][j]+BJ[i][j+1]+A[i-1][j]+A[i][j -1]-4*BJ[i][j]);//迭代公式for(i=1;i<4;i++){for(j=1;j<4;j++)if(fabs(A[i][j]-BJ[i][j])<1e-5)pp++;}N++;}while(pp<=9);pp=0;for(i=0;i<3;i++)w[M][i+1]=A[1][i+1];for(i=3;i<6;i++)w[M][i+1]=A[2][i-2];for(i=6;i<9;i++)w[M][i+1]=A[3][i-5];s[M]=N;M++;int P,Q;for(P=0;P<4;P++)for(Q=0;Q<5;Q++)A[P][Q]=0;for(P=0;P<5;P++)A[4][P]=100;}int min=s[0];int p,q;cout<<"输出所有的加速因子的迭代次数:"<<'\n'; for(q=1;q<100;q++){// cout<<s[q]<<setw(6);// if(q%12==0)// cout<<'\n';if(min>s[q]){min=s[q];p=q;}}cout<<endl;if(min==s[0])p=0;cout<<"最佳加速因子a=";cout<<d[p]<<'\n';cout<<"迭代次数为："<<min<<'\n';cout<<"最佳收敛因子对应的各个格内点的电位为："<<'\n';for( i=1;i<10;i++){cout<<w[p][i]<<'\t';if(i%3==0)cout<<'\n';}cout<<'\n';}四、程序框图迭代解程序框图五、结果分析迭代收敛的速度与α的关系最佳收敛因子的经验公式：)sin(p120πα+=（正方形场域、正方形网格）220q 1p 122+-=πα（矩形场域、正方形网格） 程序执行结果如下。

南京理工大学课程考核论文课程名称：高等数值分析论文题目：有限差分法求解偏微分方程姓名：罗晨学号：成绩：有限差分法求解偏微分方程一、主要内容1.有限差分法求解偏微分方程，偏微分方程如一般形式的一维抛物线型方程：具体求解的偏微分方程如下：2.推导五种差分格式、截断误差并分析其稳定性；3.编写MATLAB程序实现五种差分格式对偏微分方程的求解及误差分析；4.结论及完成本次实验报告的感想。

二、推导几种差分格式的过程：有限差分法（finite-differencemethods ）是一种数值方法通过有限个微分方程近似求导从而寻求微分方程的近似解。

有限差分法的基本思想是把连续的定解区域用有限个离散点构成的网格来代替；把连续定解区域上的连续变量的函数用在网格上定义的离散变量函数来近似；把原方程和定解条件中的微商用差商来近似，积分用积分和来近似，于是原微分方程和定解条件就近似地代之以代数方程组，即有限差分方程组，解此方程组就可以得到原问题在离散点上的近似解。

推导差分方程的过程中需要用到的泰勒展开公式如下：()2100000000()()()()()()()......()(())1!2!!n n n f x f x f x f x f x x x x x x x o x x n +'''=+-+-++-+-(2-1)求解区域的网格划分步长参数如下：11k k k kt t x x h τ++-=⎧⎨-=⎩(2-2) 2.1古典显格式2.1.1古典显格式的推导由泰勒展开公式将(,)u x t 对时间展开得2,(,)(,)()()(())i i k i k k k uu x t u x t t t o t t t∂=+-+-∂(2-3) 当1k t t +=时有21,112,(,)(,)()()(())(,)()()i k i k i k k k k k i k i k uu x t u x t t t o t t tuu x t o tττ+++∂=+-+-∂∂=+⋅+∂(2-4)得到对时间的一阶偏导数1,(,)(,)()=()i k i k i k u x t u x t uo t ττ+-∂+∂(2-5) 由泰勒展开公式将(,)u x t 对位置展开得223,,21(,)(,)()()()()(())2!k i k i k i i k i i u uu x t u x t x x x x o x x x x∂∂=+-+-+-∂∂(2-6)当11i i x x x x +-==和时，代入式(2-6)得2231,1,1122231,1,1121(,)(,)()()()()(())2!1(,)(,)()()()()(())2!i k i k i k i i i k i i i i i k i k i k i i i k i i i iu u u x t u x t x x x x o x x x xu u u x t u x t x x x x o x x x x ++++----⎧∂∂=+-+-+-⎪⎪∂∂⎨∂∂⎪=+-+-+-⎪∂∂⎩(2-7) 因为1k k x x h +-=，代入上式得2231,,22231,,21(,)(,)()()()2!1(,)(,)()()()2!i k i k i k i k i k i k i k i ku u u x t u x t h h o h x xu u u x t u x t h h o h x x +-⎧∂∂=+⋅+⋅+⎪⎪∂∂⎨∂∂⎪=-⋅+⋅+⎪∂∂⎩(2-8) 得到对位置的二阶偏导数2211,22(,)2(,)(,)()()i k i k i k i k u x t u x t u x t uo h x h+--+∂=+∂(2-9) 将式(2-5)、(2-9)代入一般形式的抛物线型偏微分方程得(2-10)为了方便我们可以将式(2-10)写成11122k kk k k k i i i i i i u u u u u f h ατ++-⎡⎤--+-=⎢⎥⎣⎦(2-11) ()11122k k k k k k i i i i i i u u uu u f hτατ++----+=(2-12)最后得到古典显格式的差分格式为()111(12)k k k k k i i i i i u ra u r u u f ατ++-=-+++(2-13)2r hτ=其中，古典显格式的差分格式的截断误差是2()o h τ+。

实验报告课程名称：计算方法院系：数学科学系专业班级：数应1001学号：1031110139学生姓名：姚海保指导教师：沈林开课时间：2012至2013学年第一学期一、学生撰写要求按照实验课程培养方案的要求，每门实验课程中的每一个实验项目完成后，每位参加实验的学生均须在实验教师规定的时间内独立完成一份实验报告，不得抄袭，不得缺交。

学生撰写实验报告时应严格按照本实验报告规定的内容和要求填写。

字迹工整，文字简练，数据齐全，图表规范，计算正确，分析充分、具体、定量。

二、教师评阅与装订要求1.实验报告批改要深入细致，批改过程中要发现和纠正学生实验报告中的问题，给出评语和实验报告成绩，签名并注明批改日期。

实验报告批改完成后，应采用适当的形式将学生实验报告中存在的问题及时反馈给学生。

2.实验报告成绩用百分制评定，并给出成绩评定的依据或评分标准（附于实验报告成绩登记表后）。

对迟交实验报告的学生要酌情扣分，对缺交和抄袭实验报告的学生应及时批评教育，并对该次实验报告的分数以零分处理。

对单独设课的实验课程，如学生抄袭或缺交实验报告达该课程全学期实验报告总次数三分之一以上，不得同意其参加本课程的考核。

3.各实验项目的实验报告成绩登记在实验报告成绩登记表中。

本学期实验项目全部完成后，给定实验报告综合成绩。

4.实验报告综合成绩应按课程教学大纲规定比例（一般为10-15%）计入实验课总评成绩；实验总评成绩原则上应包括考勤、实验报告、考核（操作、理论）等多方面成绩；5.实验教师每学期负责对拟存档的学生实验报告按课程、学生收齐并装订，按如下顺序装订成册：实验报告封面、实验报告成绩登记表、实验报告成绩评定依据、实验报告（按教学进度表规定的实验项目顺序排序）。

装订时统一靠左侧按“两钉三等分”原则装订。

024********-1-0.8-0.6-0.4-0.20.20.40.60.813、画出2222)sin(yxyxz++=所表示的三维曲面。

有限差分法——傅立叶稳定性分析分析差分格式稳定性的方法很多，大部份应用于线性方程，这里只介绍其中最常用的一种：傅立叶稳定性分析法。

傅立叶稳定分析法由V on Neumann 于20世纪40年代提出，所以又称为V on Neumann 稳定性分析法。

该方法的基本思想是，将解的误差作周期延拓并用傅立叶级数表示出来，然后考察每一个傅立叶级数分量的增大和衰减情况。

如果每一分量的强度（或振幅）是随时间的推移而增大的，则所讨论的差分格式是不稳定的；反之，若每一分量的振幅是随时间的推移而衰减或保持不变，则格式是稳定的。

为了进行这种分析，可以把某一分量的表达式代入到误差传播方程中，得出相邻二时间层间该分量的振幅比，通常称为放大因子。

稳定性的条件要求放大因子的绝对值（或模）小于或等于1。

当放大因子等于1时，称为中性稳定，在这种情况下任何时刻引进的误差都不会衰减或放大。

【例11.1E 】讨论逼近以下一维对流方程的FTCS 格式的稳定性：0=∂∂+∂∂xu t u α α> 0 (11.1.51)该方程的FTCS 格式为 02111=∆-+∆--++xu u t u u n i n i n i n i α (11.1.52) 将式(11.1.52)改写成易于递推计算的差分格式，有()n i n i n i n i u u u u 1112-++--=λα式中，)/(x t ∆∆λ=为网格比。

相应于上式的误差传播方程为()n i n i n i n i 1112-++--=εελαεε (11.1.53) 式中，ε是各节点上的离散量。

如果对ε在正负方向上作周期延拓，即把ε看作是以某一定值为周期的周期函数，则n ε、1+n ε可以展开为以下的傅立叶级数：⎪⎪⎩⎪⎪⎨⎧==∑∑∞-∞=++∞-∞=k kx n k n k kx n k n e C x e C x I 11I )()(εε 于是⎪⎪⎩⎪⎪⎨⎧====∑∑∞-∞=+++∞-∞=k kx n k i n n i k kx n k i n n i i i e C x e C x I 111I )()(εεεε (11.1.54) ⎪⎪⎩⎪⎪⎨⎧=±==±=∑∑∞-∞=±+++±∞-∞=±±k x x k n k i n n i k x x k n k i n n i i i e C x x e C x x )(I 1111)(I 1)()(∆∆∆εε∆εε (11.1.55) 其中，1I -=。

有限差分法一、单变量函数：用中心差分法（matlab程序见附录）计算结果如下：图1 中心差分法表1 数据对比二、一维热传导：在此取φ(x)=0，g1(t)= g2(t)=100-100*exp(-t)问题描述：已知厚度为l的无限大平板，初温0度，初始瞬间将其放于温度为100度的流体中，流体与板面间的表面传热系数为一常数。

试确定在非稳态过程中板内的温度分布。

(1)显式差分法：图3 显式差分法(2)隐式差分法：图4 隐式差分法小结：显式格式仅当时格式是稳定的。

（其中称为网格比）隐式格式从k层求k+1层时，需要求解一个阶方程组。

而且隐式格式的稳定性对网格比没有要求，即为绝对稳定的。

三、Possion方程：取f=1,R=1图5差分法图6 误差小结：观察误差曲面，其绝对误差数量级为附Matlab程序：第1题：%===========================Boundary Value Problem 1clear;clc;A=[-2.01 1 0 0 0 0 0 0 0;1 -2.01 1 0 0 0 0 0 0;0 1 -2.01 1 0 0 0 0 0;0 0 1 -2.01 1 0 0 0 0;0 0 0 1 -2.01 1 0 0 0;0 0 0 0 1 -2.01 1 0 0;0 0 0 0 0 1 -2.01 1 0;0 0 0 0 0 0 1 -2.01 1;0 0 0 0 0 0 0 1 -2.01;];c1=[0.1;0.2;0.3;0.4;0.5;0.6;0.7;0.8;0.9];C=0.01*c1-1*[0;0;0;0;0;0;0;0;1];y=A\C;x=0:0.1:1;yn=[0;y;1];ye=2*(exp(x)-exp(-x))/(exp(1)-exp(-1))-x;figure(1);plot(x,yn,'*',x,ye);legend('numerical solution','exact solution')xlabel('x','fontsize',20);ylabel('y','fontsize',20);set(gca,'fontsize',18);figure(2);err=abs(ye'-yn);plot(x,err);legend('error')xlabel('x','fontsize',20);ylabel('y','fontsize',20);set(gca,'fontsize',18);第2题：%========================Boundary Value Problem 1_Explicit %显式clear;clcl=20;%板厚h=1;%步长J=l/h;T=50;%时间tao=2.5;%步长N=T/tao;%下面组合A矩阵a=0.2;lamda=tao/(h^2);zhu=1-2*a*lamda;ce=a*lamda;a00=ones(1,J-1);a0=diag(a00);A0=zhu*a0;a01=ones(1,J-2);a1=diag(a01,1);A1=ce*a1;a2=diag(a01,-1);A2=ce*a2;A=A0+A1+A2;u(:,1)=0; %板的初始温度for i=2:N+1u(1,i)=100-100*exp(-(i-1)*tao); %边界条件u(J+1,i)=100-100*exp(-(i-1)*tao); %边界条件end% g01=u(:,1);% g02=u(:,J);for k=1:N% g01=ce*g1(1,k);% g02=ce*g2(1,k);oo=zeros(J-3,1);g(:,k)=[ce*u(1,k);oo;ce*u(J+1,k)];u(2:end-1,k+1)=A*u(2:end-1,k)+g(:,k);endt=0:h:l;x=0:tao:T;mesh(x,t,u)xlabel('t','fontsize',20);ylabel('x','fontsize',20);zlabel('T','fontsize',20);set(gca,'fontsize',18);%========================Boundary Value Problem 1_2Implicit %隐式clear;clcl=20;%板厚h=1;%步长J=l/h;T=50;%时间tao=2.5;%步长N=T/tao;%下面组合A矩阵a=0.2;lamda=tao/(h^2);zhu=1+2*a*lamda;ce=-a*lamda;a00=ones(1,J-1);a0=diag(a00);A0=zhu*a0;a01=ones(1,J-2);a1=diag(a01,1);A1=ce*a1;a2=diag(a01,-1);A2=ce*a2;A=A0+A1+A2;u(:,1)=0; %板的初始温度for i=2:N+1u(1,i)=100-100*exp(-(i-1)*tao); %边界条件u(J+1,i)=100-100*exp(-(i-1)*tao); %边界条件end% g01=u(:,1);% g02=u(:,J);for k=1:N% g01=ce*g1(1,k);% g02=ce*g2(1,k);oo=zeros(J-3,1);g(:,k)=[ce*u(1,k);oo;ce*u(J+1,k)];u(2:end-1,k+1)=inv(A)*u(2:end-1,k)-inv(A)*g(:,k);endt=0:h:l;x=0:tao:T;mesh(x,t,u)xlabel('t','fontsize',20);ylabel('x','fontsize',20);zlabel('T','fontsize',20);set(gca,'fontsize',18);第3题：%=============================used by centered difference clear;function pdemodel[pde_fig,ax]=pdeinit;pdetool('appl_cb',1);set(ax,'DataAspectRatio',[1 1 1]);set(ax,'PlotBoxAspectRatio',[1.5 1 1]);set(ax,'XLim',[-1.5 1.5]);set(ax,'YLim',[-1 1]);set(ax,'XTickMode','auto');set(ax,'YTickMode','auto');% Geometry description:pdecirc(0,0,1,'C1');set(findobj(get(pde_fig,'Children'),'Tag','PDEEval'),'String','C1')% Boundary conditions:pdetool('changemode',0)pdesetbd(4,...'dir',...1,...'1',...'0')pdesetbd(3,...'dir',...1,...'1',...'0')pdesetbd(2,...'dir',...1,...'1',...'0')pdesetbd(1,...'dir',...1,...'1',...'0')% Mesh generation:setappdata(pde_fig,'Hgrad',1.3);setappdata(pde_fig,'refinemethod','regular');pdetool('initmesh')pdetool('refine')% PDE coefficients:pdeseteq(1,...'1.0',...'0.0',...'1',...'1.0',...'0:10',...'0.0',...'0.0',...'[0 100]')setappdata(pde_fig,'currparam',...['1.0';...'0.0';...'1 ';...'1.0'])% Solve parameters:setappdata(pde_fig,'solveparam',...str2mat('0','1524','10','pdeadworst',...'0.5','longest','0','1E-4','','fixed','Inf'))% Plotflags and user data strings:setappdata(pde_fig,'plotflags',[1 1 1 1 1 1 1 1 0 0 0 1 0 0 1 0 0 1]); setappdata(pde_fig,'colstring','');setappdata(pde_fig,'arrowstring','');setappdata(pde_fig,'deformstring','');setappdata(pde_fig,'heightstring','');% Solve PDE:pdetool('solve')。

电磁场与电磁波实验报告实验项目：有限差分法一、实验目的及要求1、学习有限差分法的原理与计算步骤；2、学习用有限差分法解静电场中简单的二维静电场边值问题；3、学习用Matlab语言描述电磁场与电磁波中内容，用matlab求解问题并用图形表示出了，学习matlab语言在电磁波与电磁场中的编程思路。

二、实验内容理论学习：学习静电场中边值问题的数值法中的优先差分法的求解知识；实践学习：学习用matlab语言编写有限差分法计算二维静电场边值问题；三、实验仪器或软件Matlab7.0电脑四、实验原理有限差分法的基本思想将计算场域划分成网格，把求解场域内连续的场分布用求解网格节点上的离散数值解来代替；即用网格节点的差分方程近似代替场域内的偏微分方程来求解。

简单迭代法先对场域内的节点赋予初始值）（0,j i Φ ，这里上标（0）表示第0次近似值，即初始值。

然后再按照： ][41k 1,k ,1k 1,k ,11k ,）（）（）（）（）（++--+Φ+Φ+Φ+Φ=Φj i j i j i j i j i进行反复迭代。

若当第N 次迭代结束后，所有内节点相邻两次迭代值之间的绝对误差小于事先给定的精度，则迭代停止。

W MAX N ji N j i 〈Φ-Φ-）（）（1,, 注意：初始值的赋予是任意的；赋予初始值后，请按“从左到右、从下到上”的固定顺序依次计算各节点值； 当所有节点都算完一遍后，再用它们的新值代替旧值，即完成一次迭代。

五、实验步骤复习理论知识；编写matlab 程序；六、结果分析与问题讨论1、程序：clearX=[0,0,0,0,0;0,25,25,25,0;0,50,50,50,0;0,75,75,75,0;100,100,100,100,100]Pot=[0,0];for i=2:4for j=2:4PotX(i-1,j-1)=(X(i-1,j)+X(i,j-1)+X(i+1,j)+X(i,j+1))/4 Pot(1)=abs(PotX(i-1,j-1)-X(i,j));Pot(2)=max(Pot)endendX(2:4,2:4)=PotXnum=1;while(max(1000.*Pot)>1)Pot(2)=0; (,1,2,......) (0,1,2,......)i j k ==for i=2:4for j=2:4PotX(i-1,j-1)=(X(i-1,j)+X(i,j-1)+X(i+1,j)+X(i,j+1))/4 Pot(1)=abs(PotX(i-1,j-1)-X(i,j));Pot(2)=max(Pot)endendX(2:4,2:4)=PotXnum=num+1endsurf([0:4],[0:4],X);shading interpcolorbar('horiz')title('有限差分法计算电位图');2、运行结果X =0 0 0 0 00 25 25 25 00 50 50 50 00 75 75 75 0100 100 100 100 100%%第一次迭代PotX =18.7500Pot =6.2500 6.2500PotX =7.1440 9.8230 7.144018.7515 25.0023 18.751542.8583 52.6801 42.8583Pot =1.0e-003 *0.3815 0.7629%%第28次迭代X =0 0 0 0 00 7.1440 9.8230 7.1440 00 18.7515 25.0023 18.7515 00 42.8583 52.6801 42.8583 0100.0000 100.0000 100.0000 100.0000 100.0000num =283、波形图matlab软件在使用有限差分法研究静电场边值问题中有着重要的作用，它能够快捷有效并且准确的解决边值问题，是解决计算相对复杂问题的有效工具。

实验:有限差分法求 解接地金属槽内电 位分布

有限差分法是基于差分原理的一种数值计算 法。其基本思想想是将场域离散成很多许多 小的网格，应用差分原理，将求解连续函数 的柏松方程问题转换为求解网格节点上的差 分方程组问题。
1．1 二维柏松方程的差分格式
y
b
U0
o
a
x
导体槽中静电场的边值问题的拉普拉斯方程为
1.3 超松弛迭代法
k 1 k
i, j i, j
1 k 1 k k k (ik 4 ) 1, j i , j 1 i 1, j i , j 1 i, j
4

式中（1<α<2）为加速收敛的因子，影响 着迭1 2 2 1 3 3 2 0 h( )0 h ( 2 )0 h ( 3 )0 . y 2! y 3! y
2 3 1 1 2 3 4 0 h( )0 (h) ( 2 )0 (h) ( 3 )0 . y 2! y 3! y
4
如何实现？
将网格分成M*N列，边界点正好都是网格的节点，对所有 的节点进行编号，并记录节点的坐标位置，并用一个二维数 组进行表示u[M][N]
k 1
i, j
1 k 1 k k ik 1, j i , j 1 i 1, j
i , j 1
如何实现？
4

4项相加，且忽略高次项
2 2 1 2 3 4 40 (h) ( 2 )0 ( 2 )0 40 y x
2
得到
0
1 2 3 4
4
1.2 高斯—赛德尔迭代法
i , j
k 1

实
验
报
刘伶
告
一、实验室名称：物理电子学院仿真实验室 二、实验项目名称：电磁仿真综合实践 I 级 三、实验学时：36 学时 四、实验原理：电磁仿真的时域有限差分法 五、实验目的：
加强学生的计算机综合应用能力、尤其是运用计算机分析和解决专业问 题的能力培养，使学生对独立进行科学研究有初步实践；初步掌握一种纯数 值电磁仿真方法——时域有限差分法；初步学会综合应用一种程序设计语言 进行科学与工程计算；增强科技报告写作能力，学会相关软件使用。
11
收边界后波形就传出。 (2)终端短路时
图 a: 向 Z=L 传输时的时域波形
图 b: 向 Z=0 传输时的时域波形
图 c: 两观测面的电压值，从上到下依次是左端面和右端面
分析： 图 a 和 b 知是电磁波整个空间内随时间变化波形图的两个截图。 由仿真结果可知， 电磁波随时间的流逝在空间上是先沿 Z=L 方向传播，再沿 Z=0 方向传播。由于终端是短路 边界，因此传输进来时值是 1，经端面后反射回去的波幅度与原波形等大，且反向，符合理 论预期。图 c 是两个观测面上的电压随时间的变化情况。由波形可知两端面在电磁波传到 该处时均再现高斯脉冲，到短路边界后脉冲再一次反向返回。
常用的一阶差分格式有三种：
因为上述三种差分格式中， 只有中心差分格式精度是 2 阶， 故可用中心差分代替时间和 空间的差分：
3）Yee 氏网格
为了采用数值差分对电磁场进行计算，在空域进行 Yee 氏网格划分：
由麦克斯伟方程组可以得出电场与磁场之间的 6 组偏微分方程关系式， 用中心差分格式 代替偏导，于是可得下面 6 组差分格式——Yee 差分格式：
n 1/2 Hy
;
Ezn 1 ; Exn 1 ;

西京学数学软件实验任务书课程名称数学软件实验班级数0901学号0912020107姓名李亚强微分方程组边值问题数值算法（打靶法，有限差分法）实验课题熟悉微分方程组边值问题数值算法（打靶法，有限差实验目的分法）运用Matlab/C/C++/Java/Maple/Mathematica等其中实验要求一种语言完成微分方程组边值问题数值算法（打靶法，有限差分法）实验内容成绩教师实验二十七实验报告1、实验名称：微分方程组边值问题数值算法（打靶法，有限差分法）。

2、实验目的：进一步熟悉微分方程组边值问题数值算法（打靶法，有限差分法）。

3、实验要求：运用Matlab/C/C++/Java/Maple/Mathematica 等其中一种语言完成程序设计。

4、实验原理：1．打靶法：对于线性边值问题(1)⎩⎨⎧==∈=+'+''βα)(,)(],[)()()(b y a y b a x x f y x q y x p y 假设是一个微分算子使：L ()()Ly y p x y q x y '''=++则可得到两个微分方程：，，)(1x f Ly =α=)(1a y 0)(1='a y ，， (2)⇔)()()(111x f y x q y x p y =+'+''α=)(1a y 0)(1='a y ，，02=Ly 0)(2=a y 1)(2='a y ，， (3)⇔0)()(222=+'+''y x q y x p y 0)(2=a y 1)(2='a y 方程（2），（3）是两个二阶初值问题.假设是问题1y（2）的解，是问题（3）的解，且，则线性边值问2y 2()0y b ≠题（1）的解为: 。

1122()()()()()y b y x y x y x y b β-=+2．有限差分法：基本思想是把连续的定解区域用有限个离散点构成的网格来代替，这些离散点称作网格的节点；把连续定解区域上的连续变量的函数用在网格上定义的离散变量函数来近似；把原方程和定解条件中的微商用差商来近似， 积分用积分和来近似，于是原微分方程和定解条件就近似地代之以代数方程组，即有限差分方程组 ， 解此方程组就可以得到原问题在离散点上的近似解。

2013～2014学年第一学期实验报告
班级 学号：姓名： 实验时间：11月29日成绩：
实验
项目
一维有限差分法之稳定性分析
所属
课程
微分方程数值解
实
验
目
的
了解一维有限差分法之算法构造，如Euler法，改进Euler等。
实
验
环
境
Matlab 7.0
实
验
内
容
掌握一维有限差分法之算法构造，用Euler法、Taylor法求解课本上45页实习题，并分析稳定性
0.7277 0.7188 0.7090 0.7001 0.6903 0.6815
0.6717 0.6630 0.6532 0.6445 0.6347 0.6260
0.6162 0.6074 0.5976 0.5888 0.5789
ans =
1.0000 1.0001 1.0004 1.0006 1.0011 1.0015
1.7545 1.7346 1.7141 1.6944 1.6742 1.6547
1.6347 1.6155 1.5958 1.5769 1.5575 1.5389
1.5199 1.5017 1.4830 1.4651 1.4469 1.4294
1.4115 1.3944 1.3769 1.3602 1.3431 1.3269
u=zeros(n+1,1);
u(1)=0.5;
x(1)=0.5;
formm=1:H
x(mm+1)=x(mm)+hh*f2(mm*hh,x(mm));
end
u(2)=x(H+1);
form=1:n-1
u0=u(m)+h*f2(m*h,u(m));

有 限 差 分 法流体运动的控制方程多为偏微分方程，在复杂的情况下不存在解析解。

但是对于一些简单的情况存在解析解，偏微分方程的解析解可用精确的数学表达式表示，该表达式给出了因变量在整个定义域中的连续变化状况。

有限差分法（Finite Difference Method ,FDM ）是数值计算中比较经典的方法，由于其计算格式直观且计算简便，因此被广泛地应用在计算流体力学中。

有限差分法首先将求解区域划分为差分网格，变量信息存储在网格节点上，然后将偏微分方程的导数用差商代替，代入微分方程的边界条件，推导出关于网格节点变量的代数方程组，通过求解代数方程组，获得偏微分方程的近似解。

偏微分方程被包含离散点未知量的代数方程所替代，这个代数方程能求出离散节点处的变量，这种离散方法叫做有限差分法。

2.1有 限 差 分 逼 近2.1.1 有限差分网格 由于有限差分法求解的是网格节点上的未知量值，因此首先介绍有限差分网格。

图2.1 – 1是x-y 平面上的矩形差分网格示意图。

在x 轴方向的网格间距为△x ，在y 轴方向的网格间距为△y ，网格的交点称为节点，计算变量定义在网格节点上。

称△x 和△y 为空间步长，△x 一般不等于△y ，且△x 和△y 也可以不为常数。

取各方向等距离的网格，可以大大简化数学模型推导过程，并且经常会取得更加精确的数值解。

本章作为计算流体力学入门知识，假设沿坐标轴的各个方向网格间距分别相等，但是并不要求各方向的网格间距一致。

例如假设△x 和△y 是定值，但是不要求△x 等于△y 。

在图2.1 - 1中，网格节点在x 方向用i 表示，在y 方向用j 表示。

因此，假如（i ，j ）是点P 在图2.1 – 1中的坐标，那么，点P 右边的第一个点的就可以用（i+1，j ）表示；在P 左边的第一个点的就可以用（i —1，j ）表示；点P 上边的第一个点的就可以用（i ，j+1）表示；点P 下边的第一个点的就可以用（i ，j —1）表示。

实验一 用有限差分法解静电场边值问题一、目的1．掌握有限差分法的原理与计算步骤； 2．理解并掌握求解差分方程组的超松弛迭代法，分析加速收敛因子α的作用； 3．学会用有限差分法解简单的二维静电场边值问题，并编制计算程序。

二、方法原理有限差分法是数值计算中应用得最早而又相当简单、直观的一种方法。

应用有限差分法通常所采取的步骤是：⑴ 采用一定的网格分割方式离散化场域。

⑵ 进行差分离散化处理。

用离散的、只含有限个未知数的差分方程组，来近似代替场域内具有连续变量的偏微分方程以及边界上的边界条件（也包括场域内不同媒质分界面上的衔接条件）。

⑶ 结合选定的代数方程组的解法，编制计算机程序，求解由上面所得对应于待求边值问题的差分方程组，所得解答即为该边值问题的数值解。

现在，以静电场边值问题⎪⎩⎪⎨⎧==∂∂+∂∂)2()()1(02222s f D y x Lϕϕϕ中在为例，说明有限差分法的应用。

f (s )为边界点s 的点函数，二位场域D 和边界L 示于图5.1-1中。

x图5.1-1 有限差分的网格分割1． 离散化场域应用有限差分法时，首先需从网格划分着手决定离散点的分布方式。

通常采用完全有规律的方式，这样在每个离散点上可得出相同形式的差分方程，有效地提高解题速度。

如图5.1-1所示，现采用分别与x ，y 轴平行的等距（步距为h ）网格线把场域D 分割成足够多的正方形网格。

各个正方形的顶点（也即网格线的交点）称为网格的结点。

这样，对于场域内典型的内结点0，它与周围相邻的结点1、2、3和4构成一个所谓对称的星形。

2．差分格式造好网格后，需把上述静电场边值问题中的拉普拉斯方程（1）式离散化。

设结点0上的电位值为ϕ0。

结点1、2、3和4上的电位值相应为ϕ1、ϕ2、ϕ3和ϕ4，则基于差分原理的应用，拉普拉斯方程（1）式在结点0处可近似表达为ϕ1+ϕ2+ϕ3+ϕ4-4ϕ1=0 （3）这就是规则正方形网格内某点的电位所满足的拉普拉斯方程的差分格式，或差分方程。

实学算编实验一科学计算编程科学计算编程及及有限差分法实验计算电磁学－有限差分法电子科技大学物理电子学院赖生建cem@编程准备整程序计算路整理程序计算思路分析所需的数据分析所需的计算模块及所需数据 分析模块间的关系编写数据变量编写计算数块）编写计算函数（块）按流程连接各程序块调试cem@编程语言C语言，开发平台Microsoft Visual C++ 6.0F t语言开发平台Vi l F t60Fortran语言，开发平台Visual Fortran 6.0Matlab语言，Matlab6.5前种言科算编言编性言 前两种语言是常用的科学计算编程语言，属于编译性语言，编译后的二进制运行程序计算速度很快，按一般的编程方法，是需要先定义计算所需的各计算变量，是严格数据方法是需要先定义计算所需的各计算变量是严格类型检查性的语言，大型科学计算代码一般采用这两种语言而语言要针对矩阵计算的语言属于解释言。

而MATLAB语言主要针对矩阵计算的语言，属于性的语言，不需要先定义各计算变量，不需要数据类型检查，自动使用、识别数据变量，MATLAB具有丰富的各种科学计算工具箱和易于控制的图形库，广泛应用科学计算科学计算具箱和易于控制的图形库广泛应用科学计算中。

本实验的计算程序采用C语言编写，绘图程序采用MATLAB中的图形函数。

cem@VC6.0启动VC6.0平台C C语言编程实践-启动VC平台VCcem@启动界面VC6.0启动界面VC6.0cem@创建标准Ｃ语言程序VC6.0创建标准Ｃ语言程序VC6.0cem@在代码编辑区编写Ｃ语言程序cem@编译Ｃ语言程序cem@VC6.0调试、运行Ｃ语言程序VC6.0调试、运行Ｃ语言程序cem@平台使用MATLAB平台使用MATLABcem@的程序语言编辑MATLAB的程序语言编辑MATLABcem@科学计算数据的保存编写的计算程序运行后，计算出的数据如何处理？通常把数据保存到文件中，以备使用。

有限差分法一、有限差分法的定义有限差分法（Finite Differential Method ）是基于差分原理的一种数值计算法。

其基本思想：将场域离散为许多小网格，应用差分原理，将求解连续函数ϕ的泊松方程的问题转换为求解网格节点上ϕ的差分方程组的问题。

二、有限差分法的应用例3.7.1 有一个无限长直的金属槽，截面为正方形，两侧为正方形，两侧面及底板接地，上盖板与侧面绝缘，其上的电位为ϕ=100V, 试用有限差分法计算槽内电位。

（1）用Matlab 中的有限差分法计算槽内电位；（2）对比解析法和数值法的异同点；（3）选取一点，绘制收敛曲线；（4）总的三维电位图；1、根据有限差分公式计算出电位最终近似值为1,12,13,11,22,23,21,32,33,3=7.144=9.823=7.144=18.751=25.002=18.751=42.857=52.680=42.857ϕϕϕϕϕϕϕϕϕ，，，，，，用Matlab有限差分法计算出来结果：（见附录程序一）2、解析法和数值法的异同点解析法数值法定义在分析具体问题的基础上，抽取出一个数学模型，这个数学模型能用若干个解析表达式表示出来，解决了这些表达式，问题也就得以解决。

数值法是用高性能的计算机以数值的、程序的形式解决问题，主要是指有限元法和差分法相同点都是在具体问题的基础上取一个用解析表达式表示的数学模型来解决问题；数值法是在解析法的基础上在不同尺度上进行有限元离散,离散单元尺度不同,进行有限元计算时要满足的连续性条件不同,预测结果的精确度就不同不同点解析法可以计算出精确的数值结果；可以作为近似解和数值解的检验标准；解析法过程可以观察到问题的内在和各个参数对数值结果起的作用。

但是分析过程困难又复杂使其仅能解决很少量的问题。

数值法求解过程简单，普遍性强，用户拥有的弹性大；用户不必具备高度专业化的理论知识就可以用提供的程序解决问题。

但求解结果没有解析法精确。

利用有限差分法分析电磁场边界问题在一个电磁系统中，电场和磁场的计算对于完成该系统的有效设计师极端重要的.例如,在系统中，用一种绝缘材料是导体相互隔离是,就要保证电场强度低于绝缘介质的击穿强度。

在磁力开关中，所要求的磁场强弱，应能产生足够大的力来驱动开关。

在发射系统中进行天线的有效设计时,关于天线周围介质中电磁场分布的知识显然有实质性的意义。

为了分析电磁场，我们可以从问题所涉及的数学公式入手.依据电磁系统的特性，拉普拉斯方程和泊松方程只能适合于描述静态和准静态（低频)运行条件下的情况.但是,在高频应用中,则必须在时域或频域中求解波动方程，以做到准确地预测电场和磁场,在任何情况下，满足边界条件的一个或多个偏微分方程的解，因此,计算电池系统内部和周围的电场和磁场都是必要的。

对电磁场理论而言,计算电磁场可以为其研究提供进行复杂的数值及解析运算的方法,手段和计算结果；而电磁场理论则为计算电磁场问题提供了电磁规律，数学方程，进而验证计算结果。

常用的计算电磁场边值问题的方法主要有两大类,其每一类又包含若干种方法，第一类是解析法；第二类是数值法。

对于那些具有最简单的边界条件和几何形状规则的(如矩形、圆形等）问题，可用分离变量法和镜像法求电磁场边值问题的解析解（精确解），但是在许多实际问题中往往由于边界条件过于复杂而无法求得解析解.在这种情况下，一般借助于数值法求解电磁场的数值解。

有限差分法，微分方程和积分微分方程数值解的方法。

基本思想是把连续的定解区域用有限个离散点构成的网络来代替，这些离散点称作网格的节点；把连续定解区域上的连续变量的函数用在网格上定义的离散变量函数来近似；把原方程和定解条件中的微商用差商来近似，积分用积分和来近似,于是原微分方程和定解条件就近似地代之以代数方程组,即有限差分方程组，解此方程组就可以得到原问题在离散点上的近似解。

然后再利用插值方法便可以从离散解得到定解问题在整个区域上的近似解.差分运算的基本概念:有限差分法是指用差分来近似取代微分，从而将微分方程离散成为差分方程组。

数学物理方法课程报告题目：声波有限差分法数值模拟学生姓名：xxx学号：xxx学院：地球科学与技术学院专业班级：xxxx教师：xxx2016年 4月12日声波有限差分法数值模拟Xxx（地球科学与技术学院研15级 学号：xxx ）摘要：数值模拟是最常用的正演模拟的方法。

它通过给出的结构模型和物理参数，模拟地震波的传播轨迹，了解其规律以及过程，然后通过计算来推断观测点的地震记录。

根据求解方法，地震波方程数值解法可分为有限元法、伪谱法、有限差分法。

根据本门课程的要求，并且有限差分法具有内存占用较小，精度较高等优点，本文主要采用这种方法进行模拟。

关键词：数值模拟，声波，有限差分正文1、 引言在勘探过程中，数值模拟的作用很大。

例如：1、采集上，可用于设计或者优化野外观测系统；2、处理上，可以通过数值模拟来检验是否采用了正确的反演方法。

将正演反演不断的逼近，从而使结果更加准确；3、解释上，还可以检测一下解释的资料是否正确。

而有限差分法是数值模拟最常用的方法，本文利用有限差分法，通过对声波进行正演模拟，来了解其在地下的传播规律及特点。

2、 二维各向同性介质声波方程数值模拟 使用规则网格差分对二阶方程进行求解。

具体过程：在x 方向上，关于0x 对称分布的2N 个网格节点的坐标分别为x q x N ∆-0， x q x N ∆--10，……，x q x ∆-10，x q x ∆+10，……x q x N ∆+-10，x q x N ∆+0。

其中，x ∆表示节点间的最小间距；i q 表示任意正整数。

2N 个网格节点所对应的函 数值已知，分别为()x q x f N ∆-0，()x q x f N ∆--10，……，()x q x f ∆-10， ()x q x f ∆+10……，()x q x f N ∆+-10，()x q x f N ∆+0。

利用Taylor 级数展开求解 ()x f 在点0x 处的一阶导数近似值。