第12章MATLAB在电磁场与电磁波中的应用[MATLAB大学教程][肖汉光,邹雪,宋涛].ppt

《第12章MATLAB在电磁场与电磁波中的应用[MATLAB大学教程][肖汉光,邹雪,宋涛].ppt》由会员分享，可在线阅读，更多相关《第12章MATLAB在电磁场与电磁波中的应用[MATLAB大学教程][肖汉光,邹雪,宋涛].ppt（40页珍藏版）》请在淘文阁 - 分享文档赚钱的网站上搜索。

1、绪论o 矢量分析矢量分析o电磁场的计算与仿真电磁场的计算与仿真o电磁波的计算与仿真电磁波的计算与仿真12.1.1 矢量基本运算矢量基本运算1. 矢量加减运算矢量加减运算 两个矢量的加： 两个矢量的减： MATLAB中：数组表示向量；数组直接加减表示矢量加减；()()()xxxyyyzzz+ABABABABee+e()()()xxxyyyzzzABABABABeee【例12.1】求 和 的矢量加减运算。 23xyzAeee22yz Bee 代码如下：A=1 2 3;B=0,-2,2;AB_add=A+BAB_sub=A-B 运行结果：AB_add= 1 0 5AB_sub= 1 4 112.1.

2、1 矢量基本运算矢量基本运算2. 矢量标积和矢积矢量标积和矢积 两个矢量的标积： 两个矢量的矢积： MATLAB中：计算标积：dot(A,B);计算矢积：cross(A,B);cosxxyyzzABABABA BA BsinxyzxyzxyzAAABBBA BeeeABeA B【例12.1】求矢量 , 和 的标积和矢积。 22xyAeexyBee4zCe运行结果：AB_dot= 4AB_dot2= 4AC_dot= 0AB_cross= 0 0 0AC_cross= -8 -8 0代码如下：A=2 -2 0;B=1 -1 0;C=0 0 4;AB_dot=dot(A,B)AB_dot2=sum

3、(A.*B)AC_dot=dot(A,C)AB_cross=cross(A,B)AC_cross=cross(A,C)12.1.1 矢量基本运算矢量基本运算4. 三重积三重积 标量三重积： 矢量三重积： MATLAB中：计算标量三重积： dot(A,cross(B,C);计算矢量三重积： cross(A,cross(B,C);cossinA (BC)A B CsinsinxyzxyzxyzBBBCCCAB CeeeAB CAeA B C【例例12.4】请同学们自学！请同学们自学！12.1.2 梯度、散度和旋度的计算梯度、散度和旋度的计算 MATLAB中：计算梯度： gradient();计算散

4、度和旋度： diff();【例12.5】计算三维标量场的梯度。MATLAB代码为syms x y zf=x.*exp(-x.2 - y.2- z.2)g=gradient(f,x,y,z);X,Y,Z=meshgrid(-1:.1:1,-1:.1:1,-1:.1:1);G1=subs(g(1),x y z,X,Y,Z); G2=subs(g(2),x y z,X,Y,Z);G3=subs(g(3),x y z,X,Y,Z);quiver3(X,Y,Z,G1,G2,G3)【例12.6】计算三维矢量场的散度和旋度。MATLAB代码为syms x y zF=y*z2,x2*y*z,exp(x)*2*

5、y*z;divF=diff(F(1),x)+diff(F(2),y)+diff(F(3),z);rotF=diff(F(3),y)-diff(F(2),z),diff(F(1),z)-diff(F(3),x),diff(F(2),x)- diff(F(1),y);%绘制散度场绘制散度场X,Y,Z=meshgrid(-1.2:.2:1.2,-1:.2:1,-1:.2:1);V=subs(divF,x y z,X,Y,Z); slice(X,Y,Z,V,0,0.2,-0.8)shading interp%绘制旋度场绘制旋度场figureG1=subs(rotF(1),x y z,X,Y,Z); G

6、2=subs(rotF(2),x y z,X,Y,Z);G3=subs(rotF(3),x y z,X,Y,Z);quiver3(X,Y,Z,G1,G2,G3)三维矢量场的散度场切片图三维标量场的旋度场12.1.3 场的可视化场的可视化 MATLAB的场可视化中：三维标量分布图： mesh()、surf();二维和三维场矢量分布图：quiver()、quiver3();等值线： contour();场线：streamline()、stream2();【例12.7】画出 的标量场分布。 MATLAB代码为x,y=meshgrid(-2:.1:2,-1:.05:1);z=x.*exp(-x.2 -

7、 y.2);mesh(x,y,z);%等位线和梯度场等位线和梯度场px,py=gradient(z,.2,.15);contour(x,y,z)hold onquiver(x,y,px,py),hold offaxis image二维场的三维表示二维等值线和梯度矢量场22()( , , )xyf x y zx e12.2.1 静电场的计算与仿真静电场的计算与仿真 1、电荷源与电场的关系、电荷源与电场的关系 点电荷的电场：204rqrEe【例12.8】点电荷电场计算和仿真示例。代码如下：N=0.0005;x=-N:N/40:N;y=-N:N/40:N;X,Y=meshgrid(x,y);q=10

8、(-19);eps0=8.854187817*10(-12);E=q./(4*pi*eps0*(X.2+Y.2)+10(-19);%注意分母为零带来的计算问注意分母为零带来的计算问题题surf(X,Y,E); %电场大小曲面图电场大小曲面图xlabel(x轴轴);ylabel(y轴轴);zlabel(场强大小场强大小);set(gcf,color,w)box on;light(Position,min(x),max(y),max(max(E);shading interp;lighting phong; material shinyfigurepX,pY=gradient(E,.1,.1);%

9、梯度梯度计算计算contour(X,Y,E); %等高线等高线hold onquiver(X,Y,pX,pY); %速度场速度场xlabel(x轴轴);ylabel(y轴轴);set(gcf,color,w)hold off,axis image点电荷的电场强度大小、电场矢量场和等电位分布点电荷的电场强度大小、电场矢量场和等电位分布仿真结果如下：2利用电位计算电场利用电位计算电场电偶极子的电场仿真【例12.9】电偶极子的电场计算和仿真示例。代码如下：syms x y;eps0=8.854187817*10(-12);q=10(-19);N=10(-8);a=N/10;V1=q./(4*pi*e

10、ps0*sqrt(x-a).2+y.2);V2=-q./(4*pi*eps0*sqrt(x+a).2+y.2);V=V1+V2;E=-gradient(V); X=-N:N/37:N;Y=-N:N/37:N;X,Y=meshgrid(X,Y);Vv=double(vpa(subs(V,x,y,X,Y);Ex=double(vpa(subs(E(1),x,y,X,Y);Ey=double(vpa(subs(E(2),x,y,X,Y);AE=sqrt(Ex.2+Ey.2);Ex=Ex./AE;Ey=Ey./AE;cv=linspace(min(min(Vv),max(max(Vv),501); q

11、uiver(X,Y,Ex,Ey)xlabel(x轴轴);ylabel(y轴轴);axis tightset(gcf,color,w)figuresurf(X,Y,Vv)xlabel(x轴轴); ylabel(y轴轴); zlabel(电位电位);box on;axis tightset(gcf,color,w)light(Position,min(min(X),max(max(Y),max(max(Vv);shading interp; lighting phong; material shinyfiguresurf(X,Y,AE)xlabel(x轴轴); ylabel(y轴轴); zlabe

12、l(电场强度电场强度);box on;axis tightset(gcf,color,w)light(Position,min(min(X),max(max(Y),max(max(AE);shading interp;lighting phong;material shiny【例12.9】电偶极子的电场计算和仿真示例。仿真结果如下：电偶极子的电偶极子的(a)电位分布和电位分布和(b)电场强度分布电场强度分布电场强度矢量分布和等电位线电场强度矢量分布和等电位线12.2.2 恒定磁场的计算与仿真恒定磁场的计算与仿真 1、小电流圆环产生的恒定磁场的计算、小电流圆环产生的恒定磁场的计算【例12.10】

13、小电流圆环模型的磁场分布仿真。 设小电流圆环的半径为a，电流为I，圆环平面处于xQy平面。根据电磁场理论可以推导出小电流圆环在较远空间（相对电流圆环的尺寸而言）处的矢量磁位A的公式为式中，m为小电流圆环的磁矩，其大小为IS。根据磁感应强度与矢量磁位之间的关系 ，可得20032sin( )44Iarrm rA reB=A03032203( )2cossin 4(sincossin sincos)4(cos coscos sinsin)sin 2sincoscos2sin sincos(cossin)4rxyzxyzxyzISrISrISrB reeeeeeeeeee%离散化电流圆环的计算方法离散

14、化电流圆环的计算方法Pz=linspace(-0.01,0.01,1000);R=10-3;%圆环半径为圆环半径为1cmI=1;%电流为电流为1AN=10000;eps0=8.854187817*10(-12);mu0=4*pi*10(-7);k1=(I*mu0*Pz*R)./(2*pi*(Pz.2+R2).(3/2);k2=(I*mu0*R)./(2*pi*(Pz.2+R2).(3/2);n=0:N-1;C11=sin(2*n*pi/N);C12=cos(2*n*pi/N);C2=sin(pi/N)*ones(size(C11);Bx=k1*sum(C11.*C2);By=k1*sum(C1

15、2.*C2);Bz=k2*sum(C2);subplot(3,1,1)plot(Pz,Bx);xlabel(轴线坐标轴线坐标 /m);ylabel(B_x /T);axis tightsubplot(3,1,2)plot(Pz,By);xlabel(轴线坐标轴线坐标 /m);ylabel(B_y /T)axis tightsubplot(3,1,3)plot(Pz,Bz);xlabel(轴线坐标轴线坐标 /m);ylabel(B_z /T)axis tightset(gcf,color,white)仿真代码如下：仿真代码如下：仿真结果如下：仿真结果如下：图图12.10 电流圆环轴线上磁感应强度

16、的各个分量电流圆环轴线上磁感应强度的各个分量12.2.2 恒定磁场的计算与仿真恒定磁场的计算与仿真 2、电荷在变化磁场中的运动、电荷在变化磁场中的运动【例12.11】电荷在非均匀磁场中的运动仿真。 设电荷以一定的速度（具有x轴分量和yOz平面分量）进入均匀增强的磁场中，磁场沿着x轴方向逐渐增强，在yOz平面为均匀磁场。利用MALTAB仿真电荷非均匀的特殊磁场中的运动。MATLAB代码如下：function yundongdianhe %电荷在非均匀磁场中的运动v=10; %设定带电粒子的初速度及入射角sita=pi/6; v=v*cos(sita); %计算x，y方向的初速度u=v*sin(s

17、ita); w=0; %求解名为“yy”的微分方程组t,y=ode23(yy,0:0.002:2,0,v,0,u,0,w); figure %描绘运动轨迹plot(t,y(:,1); %绘制一般二维曲线xlabel(t);ylabel(x);figureplot(t,y(:,3);xlabel(t);ylabel(y);figureplot(t,y(:,5);xlabel(t);ylabel(z);figureplot(y(:,3),y(:,5);xlabel(y);ylabel(z);figureplot3(y(:,1),y(:,3),y(:,5) %绘制一般三维曲绘制一般三维曲线图线图bo

18、x on;comet3(y(:,1),y(:,3),y(:,5) %绘制三维动态绘制三维动态轨迹轨迹xlabel(x);ylabel(y);zlabel(z); %电荷在非均匀磁场中运动的微分方程电荷在非均匀磁场中运动的微分方程function f=yy(t,y);global A; %定义全局变量定义全局变量A=100; %设定设定qB0/mf=y(2);0;y(4);A*y(6)*y(1);y(6);-A*y(4)*y(1);%写入微分方程写入微分方程xlabel(轴线坐标轴线坐标 /m);ylabel(B_z /T)axis tightset(gcf,color,white)仿真结果如下

19、：仿真结果如下：图图12.15 电荷在非均匀磁场中的三维运动轨迹电荷在非均匀磁场中的三维运动轨迹图图12.11 电荷在电荷在x轴上的运动轨迹轴上的运动轨迹 图图12.12 电荷在电荷在y轴上的运动轨迹轴上的运动轨迹图图12.13 电荷在电荷在z轴上的运动轨迹轴上的运动轨迹 图图12.14 电荷在电荷在yz平面上的运动轨迹平面上的运动轨迹12.3.1电磁波的合成计算与仿真电磁波的合成计算与仿真振动方向相同的两个波为： 则它们的合成波为： 当 时， ，但一般条件下，合振动的解析式很难求出，但可以利用计算机模拟合成波。1111122222cos()cos()mmExEtk zExEtk z12ExE

20、xEx12cos()ExAtkz【例12.12】合成波的计算机仿真。 程序流程图如图所示，先输入振幅、频率和初相位的参数值，然后在每个时间步计算一次两个波和合成波的位置，并更新画面，从而得到动态的合成波图像。代码如下：k=2*pi;Em=20*sqrt(2);w=10;z=0:0.01:10;Z=2;for i=1:1000 t=i*0.01; Ex1=Em*cos(w*t-k*z+0.1*pi) ; Ex2=Em*cos(w*t-k*z+0.3*pi); subplot(2,1,1) plot(z,Ex1,b,z,Ex2,r);set(gcf,color,w) subplot(2,1,2)

21、plot(z,Ex1+Ex2,b) axis(0 10 -70 70);set(gcf,color,w) pause(0.1)end运行代码，即可见动态结果！运行代码，即可见动态结果！12.3.2 电磁波驻波的模拟电磁波驻波的模拟驻波由正反两个方向的行波叠加而成，正向波可表示为： 反向波可表示为： 两者叠加可形成驻波：1cos2txyAT2cos2txyAT12yyy12.3.2 电磁波驻波的模拟电磁波驻波的模拟【例12.13】驻波模拟示例。 利用MATLAB进行驻波模拟时，先设定振幅周期和波长等参数，然后利用二重循环计算两个振动每个时间步中每个质点的位置，并合成，然后对每个时间步的质点位置绘

22、图，并进行循环更新得到动态的驻波波形。程序流程图如图所示。MATLAB代码如下：k=2*pi; Em=20*sqrt(2); w=1; z=0:0.01:10for i=1:1000 t=i*0.1; Ex=Em*cos(w*t+k*z) ; Ey=1.2*Em*cos(w*t-k*z) ; subplot(2,1,1) plot(z,Ex,b,z,Ey,r);set(gcf,color,w) subplot(2,1,2) plot(z,Ex+Ey,b); hold on axis(0 10 -70 70);set(gcf,color,w) pause(0.1)end驻波的仿真结果如下：驻波的

23、仿真结果如下：运行代码，即可见动态结果！运行代码，即可见动态结果！12.3.3 光光的多缝衍射模拟的多缝衍射模拟光的多缝衍射由单缝衍射和多缝之间的干涉光强叠加而成。单缝衍射在屏上的光强分布为 多缝干涉在屏上的光强分布为： 多缝衍射在屏上的光强分布为：2211(sin/),sin/IIu uua 2220(sin/sin),sin/IINvvvd 12III【例12.14】多缝衍射模拟。 利用MATLAB进行多缝衍射模拟时，先设定单缝宽度、双缝间距、缝数和光波长等参数，然后以一定角度步长计算-90度到90度之间各角度的I1和I2，再后将两者相乘得到多缝衍射的光强分布，流程图如图所示。MATLAB

24、代码如下：a=0.5;d=2;N=4;lamda=0.1;theda=-pi/2:pi/18000+0.00001:pi/2;u=pi*a*sin(theda)/lamda;v=pi*d*sin(theda)/lamda;I1=1*(sin(u).2./(u.2);I2=2*(sin(N*v).2./(sin(v).2);I=I1.*I2;subplot(3,1,1); plot(I1);title(单缝衍射单缝衍射)axis(5000 12000 0 1.2);set(gcf,color,w)subplot(3,1,2); plot(I2);title(多缝干涉多缝干涉)axis(5000

25、12000 0 40);set(gcf,color,w)subplot(3,1,3); plot(I);title(多缝衍射多缝衍射)axis(5000 12000 0 40);set(gcf,color,w)仿真结果如下：仿真结果如下：12.3.4 电磁波的极化模拟电磁波的极化模拟两个相互正交的、频率相同、振幅不同、相位相同的线极化平面波，可以合成线极化平面波： 两个相互正交的、频率相同、振幅相同、相位相差90o的线极化平面波，可以合成圆极化平面波： 两个相互正交的、频率相同、振幅不同、相位不同的线极化平面波，可以合成椭圆极化平面波：cos()cos()xxxmyyymEtkzEtkzEeE

26、ecos()cos(/ 2)xxmyymEtkzEtkzEeEecos()cos()xxxmyyymEtkzEtkzEeEe【例12.15】电磁波的极化仿真。仿真代码如下： figure(1) plot(Ex,Ey,b.)%绘制当前点绘制当前点 hold on; %图形保持图形保持 xlabel(x轴轴); ylabel(y轴轴); axis(-30,30,-30,30,0,10); %范围范围 set(gcf,color,w)%设置背景设置背景 pause%暂停暂停 grid on;end k=2*pi;%传播常数传播常数Em=20;%电场有效值电场有效值Exm=20*sqrt(2);%电场

27、电场x分量最大值分量最大值Eym=10*sqrt(2);%电场电场y分量最大值分量最大值w=10;%角频率角频率z=0;%空间位置空间位置for i=1:1000 t=i*0.01;%时间时间 Ex=Exm*cos(w*t-k*z) ; Ey=Eym*cos(w*t-k*z) ; axis(5000 12000 0 40); set(gcf,color,w)仿真结果如下：仿真结果如下： 图12.22 合成线极化不同时刻的轨迹 图12.23 合成圆极化不同时刻的轨迹 图12.24 合成椭圆极化不同时刻的轨迹 图12.25 合成圆极化空间各点的电场在不同时刻的轨迹12.3.5 电磁电磁波传播的模拟

28、波传播的模拟1. 真空和导体中传播的电磁波真空和导体中传播的电磁波 TEM电磁波传播时，传播方向与电场和磁场相互垂直，即满足如下公式： 式中， S为能流密度矢量，代表波的传播方向。在真空中，磁场和电场同相，振幅保持不变，不发生衰减。在导体中，由于电导率不为零，传播常数为复数，导致磁场和电场不同相且振幅不断衰减。 SEH【例12.16】模拟真空和导体中电磁波的传播。仿真代码如下：figure(1) plot3(x,Ey,Zo1,b); hold on; plot3(x,Zo1,Hz,b); grid on; axis(0,3,-30,30,-30,30) xlabel(x轴轴); ylabel(

29、电场电场); zlabel(磁场磁场); set(gcf,color,w) pause(0.01) hold off;endclear;clck=2*pi;w=10;Exm=20*sqrt(2);Hym=15*sqrt(2);x=0:0.01:3;Zo1=zeros(size(x);for i=1:1000 t=i*0.01; Ey=Exm*cos(w*t-k*x) ; %真空中真空中 Hz=Hym*cos(w*t-k*x) ; %真空中真空中 %Ey=Exm*cos(w*t-k*x).*exp(-0.5*x);%导体中导体中 %Hz=Hym*cos(w*t-k*x-pi/4).*exp(-0

30、.5*x) ; %导体中导体中仿真结果如下：仿真结果如下： 图12.26 电磁波在真空中传播时不同时刻的电场和磁场 图12.27 电磁波在导体中传播时不同时刻的电场和磁场 12.3.5 电磁电磁波传播的模拟波传播的模拟2. 波导波导系统中传播的电磁波系统中传播的电磁波 矩形波导中TM模的电场强度E、磁场强度H 的场分量的表达式为：j02j02j0j02jcossinejsincosesincosejsincosejzzzzk zzxck zzyck zzk zxcxk EmmmExykaabk EnmnExykbabmnEExyabEnmnHxykbabH j02j02sincosejcoss

31、ine0zzk zck zyczEnmnxykbabEmmnHxykaabH 【例12.17】矩形波导仿真。仿真代码如下：z=0:lamda/20:2*lamda; xs=0:a/20:a; %上板上板xs,z=meshgrid(xs,z);ys=b*ones(size(xs);y2=0*ones(size(xs);z=0:lamda/20:2*lamda; %侧板侧板yc=0:b/20:b;yc,z=meshgrid(yc,z);xc=0*ones(size(yc); for t=0:T/100:5*T %上板上板 y=ys;x=xs; %下板下板 %y=y2;Ex=real(Ce*npib

32、*cos(mpia*x).*sin(npib*y).*exp(j*(w*t-kz*z); Ey=real(Ce*mpia*sin(mpia*x).*cos(npib*y).*exp(j*(w*t-kz*z); Ez=0*Ex; clear; clc;a=22.86*10(-3); b=10.16*10(-3);f=9.84*109; lamda=3*108/f;T=1/f; w=2*pi*f;m=1; n=0;eps0=8.854187817*10(-12);mu=4*pi*10(-7);kc2=(m*pi/a)2+(n*pi/b)2;kz=w*sqrt(mu*eps0);E0=106;H0=

33、E0/(sqrt(mu/eps0);Ch=j*(kz*H0/kc2);Ce=j*(w*mu*H0/kc2);mpia=m*pi/a; npib=n*pi/b;Hx=real(Ch*mpia*sin(mpia*x).*cos(npib*y).*exp(j*(w*t-kz*z); Hy=real(Ch*npib*cos(mpia*x).*sin(npib*y).*exp(j*(w*t-kz*z); Hz=real(H0*cos(mpia*x).*cos(npib*y).*exp(j*(w*t-kz*z); quiver3(z,x,y,Ez,Ex,Ey,0.5,b); hold on; mesh(z

34、,x,y) quiver3(z,x,y,Hz,Hx,Hy,0.5,m); view(-10,50) y=yc;x=xc; %侧板 Ex=real(Ce*npib*cos(mpia*x).*sin(npib*y).*exp(j*(w*t-kz*z); Ey=real(Ce*mpia*sin(mpia*x).*cos(npib*y).*exp(j*(w*t-kz*z); Ez=0*Ex; Hx=real(Ch*mpia*sin(mpia*x).*cos(npib*y).*exp(j*(w*t-kz*z); Hy=real(Ch*npib*cos(mpia*x).*sin(npib*y).*exp(j*(w*t-kz*z); Hz=real(H0*cos(mpia*x).*cos(npib*y).*exp(j*(w*t-kz*z); quiver3(z,x,y,Ez,Ex,Ey,0.5,b); mesh(z,x,y) quiver3(z,x,y,Hz,Hx,Hy,0.5,m); xlabel(Z); ylabel(X); zlabel(Y); colormap(1,1,1); axis tight set(gcf,color,w); hold off pause(0.05);end仿真结果如下：仿真结果如下： 图12.29 不同时刻矩形波导TE10模的电磁场分布图