MATLAB语言课程论文 基于MATLAB的电磁场数值图像分析

1、基于matlab的电磁场数值分析应用 摘要 matlab使用计算机进行电磁场数值分析已成为电磁场的工程开发、科研和教学的重要手段。编程实现从电磁场微分方程到有限元求解全过程需要很好的理论基础和编程技巧，难度较高。该文介绍了电磁场数值分析的基本理论并通过几个实例介绍了使用matlab 实现电磁场偏微分方程的有限元解法。 实验结果表明这一方法具有操作简单明了!运算速度快，计算误差可控制等优点关键词 电磁场数值分析 matlab 麦克斯韦方程1、 问题的提出电磁学是物理学的一个分支，是研究电场和电磁的相互作用现象。电磁学从原来互相独立的两门科学(电学、磁学)发展成为物理学中一个完整的分支学科，主要是

2、基于电流的磁效应和变化的磁场的电效应的发现。这两个实验现象，加上麦克斯韦关于变化电场产生磁场的假设，奠定了电磁学的整个理论体系，发展了对现代文明起重大影响的电工和电子技术。针对电磁场学习理论性强、概念抽象等特点，利用matlab强大的数值计算和图形技术，通过具体实例进行仿真，绘制相应的图形，使其形象化，便于对其的理解和掌握。将matlab引入电磁学中，利用其可视化功能对电磁学实验现象进行计算机模拟，可以提高学习效率于学习积极性，使学习效果明显。通过matlab软件工具，对点电荷电场、线电荷产生的电位、平面上n个电荷之间的库仑引力、仿真电荷在变化磁场中的运动等问题分别给出了直观形象的的仿真图和数

3、值分析，形实现了可视化学习，丰富了学习内容，提高了对电磁场理论知识的兴趣。从而更好地解决电磁场中数值分析的问题。二、电磁场数值解法 麦克斯韦方程组是电磁场理论的基础，也是电磁场数值分析的出发点。它的微分形式方程： （1）式中磁场强度电通密度电场强度磁感应强度。电磁场中各种场量之间的关系由媒质的特性确定。在各向同性媒质中，由下列结构方程组确定 (2)为获得电磁场问题的唯一解!除上述方程组之外尚需给出定解条件，对静态场和稳态场只需加边界条件，对时变场还需另加初始条件。边界条件包括：（1）第一类边界条件是给定边界上的值，其中是边界点的函数或常数；（2）第二类边界条件给定边界上法向导数的值；（3）第三

4、类边界条件给定边界值与法向导数的线性组合 (3)根据麦克斯韦方程组和结构方程组!在静电场中，以电位 为求解对象，在各向同性介质中，电位满足泊松方程： (4)在恒定磁场中，取矢量磁位为求解对象，令有 (5)考虑电磁波动方程： (6)在正弦稳态条件下,上式可分别导出亥姆霍兹方程: (7)如上述几个例子,对于不同的电磁场实际应用问题求解可以得到对应的电磁场偏微分方程,直接用解析法求解这些方程组往往会遇到很多困难甚至无法求解,电磁场数值分析方法已成为求解电磁场问题的重要方法。数值分析方法将原来连续的场域离散化求解， 再用离散点的结果近似逼近连续场域的解&常用的电磁场数值分析方法包括有限差分法、边界元法

5、和有限元法。 有限差分法是以差分原理为基础的一种数值计算方法，即用差分方程代替偏微分方程，把要求解的边值问题转化为一组相应的差分问题，将求解区域划分，求解差分方程组从而得出各网格单元的场值。这种方法的特点是方法简单，网格划分容易，但对不规则边界处理不便，网格划分缺乏灵活性。边界元法以麦克斯韦积分方程为基础，它采用分步积分如格林定理等在一定条件下把该积分方程转化为关于边界的积分方程，并据此进行离散获得对应的代数方程!求出场域中变量的数值。 它的特点是将数值法与解析法相结合， 在数学上起到降维的作用，减少了计算量，但对非线性情况失去了高精度特点，有局限性。有限元法由于单元定义灵活，处理边界条件容易

6、，具有正定对称系数矩阵而占据主导地位&有限元法是根据变分原理和离散化而取得近似解的方法。是先从偏微分方程边值问题出发，找出一个能量泛函的积分式，并令其在满足第一类边界条件的前提下取极值，即构成条件变分问题& 例如与上式对应的二维泊松场第三类边界条件下的泛函极值问题为： (8)与上式对应的亥姆霍兹方程泛函为： (9) 这个条件变分问题是和偏微分方程边值问题等价的。有限元法便是以条件变分问题为对象来求解电磁场问题。在求解过程中，将场的求解区域剖分成有限个单元，因此在单元中构造出插值函数，将插值函数代入能量泛函的积分式，再把泛函离散化成多元函数。通过多元函数极值求极值方法得到一个代数方程组。最后由此

7、方程组求解得到数值解。下面举几例论证： 例证1设圆线圈的中心为o，半径为r，放置于y-z平面，线圈通过的电流为i0，如右图所示。用毕奥萨伐尔定律计算载流圆线圈在z=0处x-y平面上的磁场分布。解题分析根据毕奥萨伐尔定律， (10)线圈上任一点处的电流元在x-y平面上一点p产生的元磁场为db。在编制程序时，将电流环分为n段，每一小段视为一电流元，然后求出每一电流元在观察点处的磁场分量，求出总磁场，最后叠加。 matlab程序如下： clear all %清除 r=input(请输入圆环半径,r=); %定义输入变量 i0=input(请输入电流,i0=); %定义输出变量mu0=4*pi*1e-

8、7; c0=mu0/(4*pi); %归并常数n=20; %电流环分段数x=linspace(-3,3,n); y=x; %确定观测点范围theta0=linspace(0,2*pi,n+1); %环的圆周角分段 theta1=theta0(1:n); y1=r*cos(theta1); z1=r*sin(theta1); %环各段矢量的起始坐标y1,z1theta2=theta0(2:n+1); y2=r*cos(theta2); z2=r*sin(theta2); %环各段矢量的终点坐标y2,z2xc=0; yc=(y2+y1)./2; zc=(z2+z1)./2; %计算环各段矢量中点的

9、三个坐标分量xc,yc,zcdlx=0;dly=y2-y1;dlz=z2-z1; %计算环各段矢量dl的三个长度分量，其中x1=x2=0。 ngx=n; ngy=ngx; %网格线数for i=1:ngy %循环计算各网点上的b（x,y）值 for j=1:ngxrx=x(j)-xc; ry=y(i)-yc; rz=0-zc; %计算径矢r的3个长度分量，r在z=0的面上。 r3=sqrt(rx.2+ry.2+rz.2).3; %计算r3 dlxr_x=dly.*rz-dlz.*ry; %计算叉乘dlr的x和y分量，z分量为0 dlxr_y=dlz.*rx-dlx.*rz; bx(i,j)=s

10、um(c0*i0.*dlxr_x./r3); %把环各段产生的磁场分量累加by(i,j)=sum(c0*i0.*dlxr_y./r3);b=(bx.2+by.2).0.5; %计算b的大小endendsubplot(1,2,1), quiver(x,y,bx,by), %画矢量场图hold on plot(0,1,ro,0,-1,bo),xlabel(x),ylabel(y), %修饰图形，标注坐标轴axis(-3,3,-3,3), subplot(1,2,2)mesh(x,y,b);axis(-3,3,-3,3,0,1e-4) %画磁场大小分布图xlabel(x),ylabel(y),zla

11、bel(b)运行该程序，例如，在命令窗中的r和i0的提示后分别键入1和100，运行结果如图1所示图1 载流圆线圈的磁场分布度及位置曲线例证2一对相同的圆形线圈，彼此平行而共轴。设两线圈内的电流都是i，且回绕方向一致，线圈的半径为r，二者的间距为a（当a=r时，称为亥姆霍兹线圈），求轴线附近的磁场分布。解题分析 本题是把观测区域取在两线圈之间的小范围内。线圈b生成的左边的磁场等于线圈a的左边磁场。因 此，a、b两线圈在中间部分的合成磁场等于a线圈的右磁场与其左磁场平移r后的和。matlab程序如下：clear all %清除 mu0=4*pi*1e-7;c0=mu0/(4*pi); %归并常数i

12、0=5.0;r=1; ngx=21;ngy=21; %设定网格线数x=linspace(-1,1,ngx); %确定观测点范围y=linspace(-1,1,ngy); n=20; %电流环分段数 theta0=linspace(0,2*pi,n+1); %环的圆周角分段 theta1=theta0(1:n);y1=r*cos(theta1); z1=r*sin(theta1); %环各段矢量的起始坐标y1,z1theta2=theta0(2:n+1);y2=r*cos(theta2); z2=r*sin(theta2); %环各段矢量的终点坐标y2,z2dlx=0;dly=y2-y1;dlz

13、=z2-z1; %计算环各段矢量dl的三个长度分量，其中x1=x2=0。xc=0; yc=(y2+y1)/2; zc=(z2+z1)/2;for i=1:ngy for j=1:ngx rx=x(j)-xc; ry=y(i)-yc; rz=0-zc; %计算径矢r的3个长度分量,r在z=0的面上。 r3=sqrt(rx.2+ry.2+rz.2).3; dlxr_x=dly.*rz-dlz.*ry; %计算叉乘dlr的x和y分量，z分量为0 dlxr_y=dly.*rx-dlx.*rz; bx(i,j)=sum(c0*i0*dlxr_x./r3); %把环各段产生的磁场分量累加 by(i,j)=

14、sum(c0*i0*dlxr_y./r3);endendbax=bx(:,11:21)+bx(:,1:11); %把x0区域bay=by(:,11:21)+by(:,1:11);subplot(1,2,1),mesh(x(11:21),y,bax);xlabel(x);ylabel(y);zlabel(b);subplot(1,2,2),quiver(x,y,bx,by,1.5), axis(square),axis(-1,1,-1,1),xlabel(x);ylabel(y); 运行结果如图2所示。可以看出，在轴线附近磁场大小均匀且沿x方向。图2 亥姆霍兹线圈轴线附近bx在x-y平面上的分布

15、及矢量场例证3设带电粒子质量为m，带电量为，电场强度e沿方向，磁感应强度b沿qyz方向. 则带电粒子在均匀电磁场中的运动微分方程为 (11)令， 则上面微分方程可化作： (12)选择和为参量，就可以分别研究0e，0=b和，等情况. 编写matlab程序如下：clear %清除syms w x y z t b e m q; %定义变量e=input(e=);b=input(b=); %输入e和b值x,y,z=dsolve(d2x=q*b/m*dy,d2y=q*e/m-q*b/m*dx,d2z=0,x(0)=0,y(0)=0,z(0)=0,dx(0)=0.01,dy(0)=6,dz(0)=0.01

16、) ; %初始条件取x(0)=y(0)=z(0)=0,dx(0)=0.01,dy(0)=6,dz(0)=0.01q=1.6e-2; m=0.02; %赋值x=subs(x y z); x=x(1),y=x(2),z=x(3), ezplot3(x(1),x(2),x(3) %绘图函数调用运行上述程序，例如，取e=4, b=8可得下列特解并给出图运行结果如图3所示。 研究时可以采用不同的初始条件和不同的参量观察不同的现象。运算特解如下：x =(t*exp(1)/8 - (100*i*exp(1) - 8*i + 4800)/(10240*exp(32*i*t)/5) - (exp(32*i*t)

17、/5)*(8*i - 100*i*exp(1) + 4800)/10240 + 15/16y =(5*exp(1)/256 + (i*(100*i*exp(1) - 8*i + 4800)/(10240*exp(32*i*t)/5) - (i*exp(32*i*t)/5)*(8*i - 100*i*exp(1) + 4800)/10240 - 1/640z =t/100图3现有e=4, b=8参数运行结果例证4如右图所示，求垂直于无限长载流直导线的平面内磁感应强度的分布。解题分析设场点p的位置为，电流元位置为，电流元矢量为。由此，场点p相对于电流元的位置矢量为 利用行列式计算idlr ，可写为

18、 （13）也可利用matlab中的det 命令函数来求该行列式，matlab程序如下：syms dx dy dz x0 x y0 y z0 z; %定义变量dl=dx,dy,dz; %定义行列式r=x0-x,y0-y, z0-z;d1cr=cross(dl,r) %求dlr的积运行结果为d1cr = dy*(z0-z)-dz*(y0-y), dz*(x0-x)-dx*(z0-z), dx*(y0-y)-dy*(x0-x)即又，r的大小为 设载流导体通过坐标原点垂直于 x-y平面放置，电流元 idl沿z轴正向，场点p位于x-y平面上。对本题目而言，dx=dy=0,x=y=0, z0=0, 矢量叉

19、乘积为，r的大小为 由毕奥萨伐尔定律 （14）有； ； （15） matlab程序 1、 用符号运算求b的表达式syms c0 i z x y r r0; %定义变量bx=c0.*i.*int(-y./(x.2+y.2+z.2).(3/2),z,-inf,inf)by=c0.*i.*int(x./(x.2+y.2+z.2).(3/2),z,-inf,inf)b=(bx.2+by.2).0.5 运行结果：bx =-2*c0*i*y/(x2+y2)(3/2)/(1/(x2+y2)(1/2)by =2*c0*i*x/(x2+y2)(3/2)/(1/(x2+y2)(1/2)b =(4*c02*i2*y

20、2/(x2+y2)2+4*c02*i2*x2/(x2+y2)2)(1/2) 即 ； （16）2、绘制磁场大小分布图和矢量场图x=-0.5:0.05:0.5;y=x; %赋值i=input(请输入电流i=); %设置输入mu0=4*pi*1e-7; c0=mu0/(4*pi); %设置函数x,y=meshgrid(x,y);bx =-2.*c0.*i.*y./(x.2+y.2).(3./2.)./(1./(x.2+y.2).(1./2);by =2.*c0.*i.*x./(x.2+y.2).(3./2)./(1./(x.2+y.2).(1./2);b=(4.*c0.2.*i.2.*y.2./(x

21、.2+y.2).2+4.*c0.2.*i.2.*x.2./(x.2+y.2).2).(1./2); subplot(1,2,1) %选择12区域中的1号区quiver(x,y,bx,by,2), axis(-0.5,0.5,-0.5,0.5), axis(square),subplot(1,2,2) %选择12区域中的2号区mesh(x,y,b) %绘图运行该程序，在命令窗中的提示后键入i0值 (例如，取i0100a)，便得到图4所示图形。图4长载流导线的磁场在x-y平面上的分布三、结论从以上利用matlab语言对几种电磁场模型的分析我们不难的出以下结论：电磁场与电磁波理论作为电子信息类专业的一门重要