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

主要内容主要内容第第1 1章章 MATLABMATLAB概述及系统环境概述及系统环境第第2 2章章 MATLABMATLAB矩阵及其运算矩阵及其运算第第3 3章章 MATLABMATLAB数据可视化数据可视化第第4 4章章 MATLABMATLAB程序设计程序设计第第5 5章章 MATLABMATLAB数值计算数值计算第第6 6章章 MATLABMATLAB符号运算符号运算第第7 7章章 MATLAB GUIMATLAB GUI设计设计第第8 8章章 MATLAB SimulinkMATLAB Simulink仿真仿真第第9 9章章 MATLABMATLAB在电路仿真中的应用在电路仿真中的应用第第1010章章 MATLABMATLAB在数字信号处理中的应用在数字信号处理中的应用绪论第第1111章章 MATLABMATLAB在数字图像处理中的应用在数字图像处理中的应用第第1212章章 MATLABMATLAB在电磁场与电磁波中的应用在电磁场与电磁波中的应用第十二章第十二章 MATLABMATLAB在电磁场与电磁波中的应用在电磁场与电磁波中的应用o 矢量分析矢量分析o电磁场的计算与仿真电磁场的计算与仿真o电磁波的计算与仿真电磁波的计算与仿真12.1 矢量分析12.1.1 矢量基本运算矢量基本运算1.矢量加减运算矢量加减运算 两个矢量的加：两个矢量的减：MATLAB中：数组表示向量；数组直接加减表示矢量加减；12.1 矢量分析【例12.1】求 和 的矢量加减运算。代码如下：A=123;B=0,-2,2;AB_add=A+BAB_sub=A-B运行结果：AB_add=105AB_sub=14112.1 矢量分析12.1.1 矢量基本运算矢量基本运算2.矢量标积和矢积矢量标积和矢积 两个矢量的标积：两个矢量的矢积：MATLAB中：计算标积：dot(A,B);计算矢积：cross(A,B);12.1 矢量分析【例12.1】求矢量 ,和 的标积和矢积。运行结果：AB_dot=4AB_dot2=4AC_dot=0AB_cross=000AC_cross=-8-80代码如下：A=2-20;B=1-10;C=004;AB_dot=dot(A,B)AB_dot2=sum(A.*B)AC_dot=dot(A,C)AB_cross=cross(A,B)AC_cross=cross(A,C)12.1 矢量分析12.1.1 矢量基本运算矢量基本运算4.三重积三重积 标量三重积：矢量三重积：MATLAB中：计算标量三重积：dot(A,cross(B,C);计算矢量三重积：cross(A,cross(B,C);【例例12.4】请同学们自学！请同学们自学！12.1 矢量分析12.1.2 梯度、散度和旋度的计算梯度、散度和旋度的计算 MATLAB中：计算梯度：gradient();计算散度和旋度：diff();【例12.5】计算三维标量场的梯度。MATLAB代码为symsxyzf=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),xyz,X,Y,Z);G2=subs(g(2),xyz,X,Y,Z);G3=subs(g(3),xyz,X,Y,Z);quiver3(X,Y,Z,G1,G2,G3)12.1 矢量分析【例12.6】计算三维矢量场的散度和旋度。MATLAB代码为symsxyzF=y*z2,x2*y*z,exp(x)*2*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,xyz,X,Y,Z);slice(X,Y,Z,V,0,0.2,-0.8)shadinginterp%绘制旋度场绘制旋度场figureG1=subs(rotF(1),xyz,X,Y,Z);G2=subs(rotF(2),xyz,X,Y,Z);G3=subs(rotF(3),xyz,X,Y,Z);quiver3(X,Y,Z,G1,G2,G3)三维矢量场的散度场切片图三维标量场的旋度场12.1 矢量分析12.1.3 场的可视化场的可视化 MATLAB的场可视化中：三维标量分布图：mesh()、surf();二维和三维场矢量分布图：quiver()、quiver3();等值线：contour();场线：streamline()、stream2();12.1 矢量分析【例12.7】画出的标量场分布。MATLAB代码为x,y=meshgrid(-2:.1:2,-1:.05:1);z=x.*exp(-x.2-y.2);mesh(x,y,z);%等位线和梯度场等位线和梯度场px,py=gradient(z,.2,.15);contour(x,y,z)holdonquiver(x,y,px,py),holdoffaxisimage二维场的三维表示二维等值线和梯度矢量场12.2 电磁场的计算与仿真12.2.1 静电场的计算与仿真静电场的计算与仿真 1、电荷源与电场的关系、电荷源与电场的关系 点电荷的电场：【例12.8】点电荷电场计算和仿真示例。代码如下：N=0.0005;x=-N:N/40:N;y=-N:N/40:N;X,Y=meshgrid(x,y);q=10(-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)boxon;light(Position,min(x),max(y),max(max(E);shadinginterp;lightingphong;materialshinyfigurepX,pY=gradient(E,.1,.1);%梯度梯度计算计算contour(X,Y,E);%等高线等高线holdonquiver(X,Y,pX,pY);%速度场速度场xlabel(x轴轴);ylabel(y轴轴);set(gcf,color,w)holdoff,axisimage12.2 电磁场的计算与仿真点电荷的电场强度大小、电场矢量场和等电位分布点电荷的电场强度大小、电场矢量场和等电位分布仿真结果如下：12.2 电磁场的计算与仿真2利用电位计算电场利用电位计算电场电偶极子的电场仿真【例12.9】电偶极子的电场计算和仿真示例。代码如下：symsxy;eps0=8.854187817*10(-12);q=10(-19);N=10(-8);a=N/10;V1=q./(4*pi*eps0*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);quiver(X,Y,Ex,Ey)xlabel(x轴轴);ylabel(y轴轴);axistightset(gcf,color,w)figuresurf(X,Y,Vv)xlabel(x轴轴);ylabel(y轴轴);zlabel(电位电位);boxon;axistightset(gcf,color,w)light(Position,min(min(X),max(max(Y),max(max(Vv);shadinginterp;lightingphong;materialshinyfiguresurf(X,Y,AE)xlabel(x轴轴);ylabel(y轴轴);zlabel(电场强度电场强度);boxon;axistightset(gcf,color,w)light(Position,min(min(X),max(max(Y),max(max(AE);shadinginterp;lightingphong;materialshiny12.2 电磁场的计算与仿真【例12.9】电偶极子的电场计算和仿真示例。仿真结果如下：电偶极子的电偶极子的(a)电位分布和电位分布和(b)电场强度分布电场强度分布电场强度矢量分布和等电位线电场强度矢量分布和等电位线12.2 电磁场的计算与仿真12.2.2 恒定磁场的计算与仿真恒定磁场的计算与仿真 1、小电流圆环产生的恒定磁场的计算、小电流圆环产生的恒定磁场的计算【例12.10】小电流圆环模型的磁场分布仿真。设小电流圆环的半径为a，电流为I，圆环平面处于xQy平面。根据电磁场理论可以推导出小电流圆环在较远空间（相对电流圆环的尺寸而言）处的矢量磁位A的公式为式中，m为小电流圆环的磁矩，其大小为IS。根据磁感应强度与矢量磁位之间的关系 ，可得12.2 电磁场的计算与仿真%离散化电流圆环的计算方法离散化电流圆环的计算方法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(C12.*C2);Bz=k2*sum(C2);subplot(3,1,1)plot(Pz,Bx);xlabel(轴线坐标轴线坐标/m);ylabel(B_x/T);axistightsubplot(3,1,2)plot(Pz,By);xlabel(轴线坐标轴线坐标/m);ylabel(B_y/T)axistightsubplot(3,1,3)plot(Pz,Bz);xlabel(轴线坐标轴线坐标/m);ylabel(B_z/T)axistightset(gcf,color,white)仿真代码如下：仿真代码如下：12.2 电磁场的计算与仿真仿真结果如下：仿真结果如下：图图12.10 电电流流圆环轴线圆环轴线上磁感上磁感应应强强度的各个分量度的各个分量12.2 电磁场的计算与仿真12.2.2 恒定磁场的计算与仿真恒定磁场的计算与仿真 2、电荷在变化磁场中的运动、电荷在变化磁场中的运动【例12.11】电荷在非均匀磁场中的运动仿真。设电荷以一定的速度（具有x轴分量和yOz平面分量）进入均匀增强的磁场中，磁场沿着x轴方向逐渐增强，在yOz平面为均匀磁场。利用MALTAB仿真电荷非均匀的特殊磁场中的运动。MATLAB代码如下：functionyundongdianhe%电荷在非均匀磁场中的运动v=10;%设定带电粒子的初速度及入射角sita=pi/6;v=v*cos(sita);%计算x，y方向的初速度u=v*sin(sita);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);12.2 电磁场的计算与仿真figureplot(t,y(:,5);xlabel(t);ylabel(z);figureplot(y(:,3),y(:,5);xlabel(y);ylabel(z);figureplot3(y(:,1),y(:,3),y(:,5)%绘制一般三维曲绘制一般三维曲线图线图boxon;comet3(y(:,1),y(:,3),y(:,5)%绘制三维动态绘制三维动态轨迹轨迹xlabel(x);ylabel(y);zlabel(z);%电荷在非均匀磁场中运动的微分方程电荷在非均匀磁场中运动的微分方程functionf=yy(t,y);globalA;%定义全局变量定义全局变量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)axistightset(gcf,color,white)12.2 电磁场的计算与仿真仿真结果如下：仿真结果如下：图图12.15 电电荷在非均匀磁荷在非均匀磁场场中的三中的三维维运运动轨动轨迹迹图图12.11 电电荷在荷在x轴轴上的运上的运动轨动轨迹迹 图图12.12 电电荷在荷在y轴轴上的运上的运动轨动轨迹迹图图12.13 电电荷在荷在z轴轴上的运上的运动轨动轨迹迹 图图12.14 电电荷在荷在yz平面上的运平面上的运动轨动轨迹迹12.3 电磁波的计算与仿真12.3.1电磁波的合成计算与仿真电磁波的合成计算与仿真振动方向相同的两个波为：则它们的合成波为：当 时，但一般条件下，合振动的解析式很难求出，但可以利用计算机模拟合成波。12.3 电磁波的计算与仿真【例12.12】合成波的计算机仿真。程序流程图如图所示，先输入振幅、频率和初相位的参数值，然后在每个时间步计算一次两个波和合成波的位置，并更新画面，从而得到动态的合成波图像。代码如下：k=2*pi;Em=20*sqrt(2);w=10;z=0:0.01:10;Z=2;fori=1:1000t=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)plot(z,Ex1+Ex2,b)axis(010-7070);set(gcf,color,w)pause(0.1)end运行代码，即可见动态结果！运行代码，即可见动态结果！12.3 电磁波的计算与仿真12.3.2 电磁波驻波的模拟电磁波驻波的模拟驻波由正反两个方向的行波叠加而成，正向波可表示为：反向波可表示为：两者叠加可形成驻波：12.3 电磁波的计算与仿真12.3.2 电磁波驻波的模拟电磁波驻波的模拟【例12.13】驻波模拟示例。利用MATLAB进行驻波模拟时，先设定振幅周期和波长等参数，然后利用二重循环计算两个振动每个时间步中每个质点的位置，并合成，然后对每个时间步的质点位置绘图，并进行循环更新得到动态的驻波波形。程序流程图如图所示。MATLAB代码如下：k=2*pi;Em=20*sqrt(2);w=1;z=0:0.01:10fori=1:1000t=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);holdonaxis(010-7070);set(gcf,color,w)pause(0.1)end12.3 电磁波的计算与仿真驻波的仿真结果如下：驻波的仿真结果如下：运行代码，即可见动态结果！运行代码，即可见动态结果！12.3 电磁波的计算与仿真12.3.3 光光的多缝衍射模拟的多缝衍射模拟光的多缝衍射由单缝衍射和多缝之间的干涉光强叠加而成。单缝衍射在屏上的光强分布为 多缝干涉在屏上的光强分布为：多缝衍射在屏上的光强分布为：12.3 电磁波的计算与仿真【例12.14】多缝衍射模拟。利用MATLAB进行多缝衍射模拟时，先设定单缝宽度、双缝间距、缝数和光波长等参数，然后以一定角度步长计算-90度到90度之间各角度的I1和I2，再后将两者相乘得到多缝衍射的光强分布，流程图如图所示。MATLAB代码如下：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(50001200001.2);set(gcf,color,w)subplot(3,1,2);plot(I2);title(多缝干涉多缝干涉)axis(500012000040);set(gcf,color,w)subplot(3,1,3);plot(I);title(多缝衍射多缝衍射)axis(500012000040);set(gcf,color,w)12.3 电磁波的计算与仿真仿真结果如下：仿真结果如下：12.3 电磁波的计算与仿真12.3.4 电磁波的极化模拟电磁波的极化模拟两个相互正交的、频率相同、振幅不同、相位相同的线极化平面波，可以合成线极化平面波：两个相互正交的、频率相同、振幅相同、相位相差90o的线极化平面波，可以合成圆极化平面波：两个相互正交的、频率相同、振幅不同、相位不同的线极化平面波，可以合成椭圆极化平面波：12.3 电磁波的计算与仿真【例12.15】电磁波的极化仿真。仿真代码如下：figure(1)plot(Ex,Ey,b.)%绘制当前点绘制当前点holdon;%图形保持图形保持xlabel(x轴轴);ylabel(y轴轴);axis(-30,30,-30,30,0,10);%范围范围set(gcf,color,w)%设置背设置背景景pause%暂停暂停gridon;endk=2*pi;%传播常数传播常数Em=20;%电场有效值电场有效值Exm=20*sqrt(2);%电场电场x分量最大值分量最大值Eym=10*sqrt(2);%电场电场y分量最大值分量最大值w=10;%角频率角频率z=0;%空间位置空间位置fori=1:1000t=i*0.01;%时间时间Ex=Exm*cos(w*t-k*z);Ey=Eym*cos(w*t-k*z);axis(500012000040);set(gcf,color,w)12.3 电磁波的计算与仿真仿真结果如下：仿真结果如下：图12.22 合成线极化不同时刻的轨迹 图12.23 合成圆极化不同时刻的轨迹12.3 电磁波的计算与仿真 图12.24 合成椭圆极化不同时刻的轨迹 图12.25 合成圆极化空间各点的电场在不同时刻的轨迹12.3 电磁波的计算与仿真12.3.5 电磁电磁波传播的模拟波传播的模拟1.真空和导体中传播的电磁波真空和导体中传播的电磁波 TEM电磁波传播时，传播方向与电场和磁场相互垂直，即满足如下公式：式中，S为能流密度矢量，代表波的传播方向。在真空中，磁场和电场同相，振幅保持不变，不发生衰减。在导体中，由于电导率不为零，传播常数为复数，导致磁场和电场不同相且振幅不断衰减。12.3 电磁波的计算与仿真【例12.16】模拟真空和导体中电磁波的传播。仿真代码如下：figure(1)plot3(x,Ey,Zo1,b);holdon;plot3(x,Zo1,Hz,b);gridon;axis(0,3,-30,30,-30,30)xlabel(x轴轴);ylabel(电场电场);zlabel(磁场磁场);set(gcf,color,w)pause(0.01)holdoff;endclear;clck=2*pi;w=10;Exm=20*sqrt(2);Hym=15*sqrt(2);x=0:0.01:3;Zo1=zeros(size(x);fori=1:1000t=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.5*x);%导体中导体中12.3 电磁波的计算与仿真仿真结果如下：仿真结果如下：图12.26 电磁波在真空中传播时不同时刻的电场和磁场 图12.27 电磁波在导体中传播时不同时刻的电场和磁场 12.3 电磁波的计算与仿真12.3.5 电磁电磁波传播的模拟波传播的模拟2.波导波导系统中传播的电磁波系统中传播的电磁波 矩形波导中TM模的电场强度E、磁场强度H 的场分量的表达式为：12.3 电磁波的计算与仿真【例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);fort=0:T/100:5*T%上板上板y=ys;x=xs;%下板下板%y=y2;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;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=E0/(sqrt(mu/eps0);Ch=j*(kz*H0/kc2);Ce=j*(w*mu*H0/kc2);mpia=m*pi/a;npib=n*pi/b;12.3 电磁波的计算与仿真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);holdon;mesh(z,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);axistightset(gcf,color,w);holdoffpause(0.05);end12.3 电磁波的计算与仿真仿真结果如下：仿真结果如下：图12.29 不同时刻矩形波导TE10模的电磁场分布图