整理平板传输线电磁仿真

1、电磁仿真综合实践NumericalSimulationofElectromagneticField作者张泳AuthorYongZhang指导老师王秉中InstructorBing-ZhongWang单位电子科大物理电子学院99777Class99777CPEofUESTC设计时间2002年5月8日一6月13日【摘要】仿真是通过建立系统模型，并且对模型进行实际研究，以代替对于实际系统的研究。Matlab5.3是一种功能强大工程数学软件。本文阐明了怎样利用Matlab5.3实现一维以及二维电磁仿真实例的建模分析、编程调试、结果分析。AbstractSimulationisanavailableme

2、ansinstudyingrealsystembyModel-building.ThispaperbuildinganumericalsimulationofelectromagneticfieldmodelbasedontheMatlab5.3,thedesignisdiscussedindetail.【关键词】电磁仿真，计算电磁学，电磁场KeyWordsElectromagneticSimulation,ComputationalElectromagnetics,ElectromagneticField【设计目的】1 .课程性质：属于非计算机专业计算机基础教育第三层次一计算机综合应用层次的

3、课程2 .任务：2.1加强学生的计算机综合应用能力、尤其是运用计算机分析和解决专业问题的能力培养，使学生对独立进行科学研究有初步实践。2.2初步掌握一种纯数值电磁仿真方法一时域有限差分法。2.3初步学会综合应用一种程序设计语言进行科学与工程计算。2.4增强科技报告的写作能力，学会相关软件使用正文第一部分：电磁仿真的意义现代科学研究的基本模式：“科学实验，理论分析，数值计算”三位一体。计算电磁学：以电磁场理论为基础，以计算机及计算机技术为工具和手段，运用计算科学所提供的各种方法，解决复杂电磁问题的一门应用科学。【电场问题求解方法】解析法一电磁场问题的传统解法。强调电磁分析和数学分析通常给出紧凑高

4、效的计算程序。程序的最终拥护只具有很少的弹性，只能改变很少的参数，结构的主要特性已被编进程序，适用于专用程序开发。纯数值法一直接以数值的、程序的形式描述电磁场问题。普适性强，用户弹性大。特定问题的边界条件、电气结构、激励等特性可以不编入基本程序，而由用户输入，更好的情况是通过图形界面输入。用户不必具备高深的电磁理论、数学及数值技术方面的知识。适用于通用商业程序开发。受到硬件限制大。FDTDTLMFEM第二部分时域有限差分法基本原理【差分原理】dff(x)f(xh)f(x)刖项差分dxxh尸工/公dff(x)f(x)f(xh)dxxh中心差分5f(x)f(xh)f(xh)dxx2h用中心差分代替

5、对空间、时间坐标的微分n1n1F-ijk)F(i2k)F(i/，k)xx误差qh)误差qh)误差qh*h)(x)2)Fn(i,j,k)x【Yee的差分网格】iin-n-F2(i,j,k)F2(i,j,k)t2(t)每个磁场分量被四个电场分量环绕，每个电场分量被四个磁场分量环绕。E和H的空间取样位置相差半个空间步长，E和H时间取样位置相差半个时间步长。Yee的差分算法考虑空间一个无缘区域，其媒质的参数不随时间变化且各项同性，则麦克斯韦旋度方程在直角坐标系中，写成分量式：Hx=1_EyEztzyHy1EzExtxzHz1ExEytyxEx1HzHytyzEy1ExHztzxEz1HyHxtzy(6

6、)【蛙跳格式】1一,v丁1.【解的稳定条件】：t【数值色散】2(y)2(z)2满足min/20min/10。第三部分均匀平行板传输线传输特性仿真在时域有限差分网格中,数值波模的传播速度将随频率改变，即有色散。这种色散由数值网格引起，而非物理上客观存在。为了减小色散，实用中通常取空间步长基本过程【平行板传输线中的主模TEM1根据平行板传输线中TEM波的特性对MAXWEL旋度方程简化得:Hyt1ExExt1Hy其中Ex(x,y,z,t)E°F(t)F(t)Hy(x,y,z,t)且F(t)F(t)0若取时域波形为高斯波形，F(t)exptt02/T2,如下图所示。zL【一维差分格式】N1/

7、2yk1/2E：1(k)E：(k)1/2tF：(k)e1(K1)K1/2JELJ：ztH；1/2(k1/2)H；1/2(k1/2)zExHyExk=Nz-1/2k=Nz【网格划分】ExHyExHyExHy4ff1-eI-e+&k=0k=1/2k=1k=3/2k=21T=0.5ns,则fmax亓=1GHz,d=0.18m,L=6m,min=0.3m,zq=0.015m,Nz400,t20z2c【时间步进过程】在初始(n=0)时令所有场两为零。将高斯脉冲激励加在输入口E：Fn1n一由FDTD方程算出n+1/2时刻所有场点的磁场H2由FDTD方程算出n+1时刻所有场点的电场En1令金属边界上

8、的切向电场分量为零，并在截断边界上使用截断边界条件记录并输出两个观测面处的电压值Vref1(nt),Vref2(nt)nn+1重复上一步骤，直到(n=nt)脉冲波全部传出计算区域。【截断边界条件】终端短路：Ex(zL,t)QE；(Nz1)0终端匹配：(吸收边界条件)Un1(0)Un(1)JtzUn1Un(0)ctzUn1(Nz)Un(Nz1)c-tzUn1(Nz1)Un(Nz)ctz【一维仿真实例】终端匹配：入射波大约经过1000步后到达边界，由于采用吸收边界条件，高斯波几乎被完全吸收，通过Ex图可以观测到吸收后只有Ex振幅小于104的波形存在，这在误差的允许范围之内。第一条、第二条分别为第一

9、、二观测面的电压。-4wX14电场强度V/M电场强度V/电压强度V0.40.3,0.20.1,0,I050100150200250300350400450Z方向网格数(0.015m)入射波_66ILLLJJJ050100150200250300350400450Z方向网格数(0.015m)被吸收后的波形_0.(01002003004005006007008009001000观测面电压时间步长2.5*1011S终端短路：入射波到达终端后由于截断边界条件采取的是短路，所以波形反转，向入射口传输，由于Nt取得足够大(Nt=2000),波形会到达入射口，所以应该在波形到达入射口之前，激励脉冲已完全进入

10、入射口后，将入射口边界条件切换为吸收边界条件。第一条第二条分别为第一、二观测面的电压。在横轴下方的是反射电压，上方为入射电压。-4o电场强度V/O-0-0M电场强度V/电压强度V_0.6_0.8-1-111hh-050100150200250300350400450Z方向网格数(0.015m)-(r1-8-10-12-14111111050100150200250300350400450Z方向网格数(0.015m)-0.1-0.-0.-0.0200400600800100012001400160018002000时间步长60-62.5*1011S被反射的波形被吸收后的波形观测面电压第四部分带挡

11、板的平行板传输线传输特性仿真【二维FDTDT程】由于不连续性的出现，沿X方向场的均匀性遭到破坏,沿Y方向场的均匀性仍然保持，此时平行板传输线中的场分量有Ex,Ez,Hy.n1Hy21n.Hy2i1,k-2xN1Ezi,k2EXi9kExni9k1zE；E；nHy12,K1n2Hy2yzE；1i,k2E；i,knHy【参数设置】:数组EX(NX,NZ+1),EZ(NX+1,NZ),HY(NX,NZ),L/2应选足够长，使得在挡板反射波返回入射口之前,边界条件。激励脉冲已完全进入入射口,并且入射口边界条件已切换为吸收【观测面电压波形】Vref1nVref2ndExz0dExz0Nref1z,ntd

12、xdNref2z,ntdxd【入射、反射、透射电压频谱】j2ftNttj2dtVt0aj2ft,edtExzExzNref1z,ntNref2z,ntNtj2fnttn1Vnte【网络的S参数】现在我们再记录下总电压Vrf：(t)、Vrf2(t)、就可以分离出入射电压Vref1(t)、反射V,；"t)、透射电压V，(t)。retotalintrtotalVref1(t)=Vref1-Vref1(t)Vref2(t)=Vref2(t)终端短路的反射系数：(f)Vre；1(f)/Vre；1(f)终端匹配的S参数：31(f)Vree1(f)/Vr：f1(f)S21(f)Vr：f2(f)/V

13、ref1(f)【终端条件】终端匹配ln1Ex1:nx,nzn1Ex1:nx,1nctzn1n,1Ex1:nx,nzEx1:nx,nzEx1:nx,nzctznctzn1n,Ex1:nx,2Ex1:nx,2Ex1:nx,1ctz终端短路E；11:nx,nz1=0【观测面电压波形】dVref1nt0ExzdVref2nt0ExzNref1乙ntdxdNref2乙ntdxdExzNref1z,ntExzNref2z,nt【入射、反射、透射电压频谱】j2ftNttj2ftNtj2fntVfVtej2dtVtej2dtt'nte0n1【二维仿真实例】终端匹配：可以通过对Ex波形的观测知道，波形到

14、达挡板时产生尖端效应，高斯波有一部分被“挤”过挡板，有一部分被挡板反射。此时应该将入射口和终端都设为吸收边界条件o如下图所示:电场强度ExV/M电场强度ExV/M(0.015m)(0.015m)(0.015m)入射波到达挡板(0.015m)(0.015m)(0.015m)(0.015m)(0.015m)Z方向网格数Z方向网格数X方向网格数电场强度ExV/MX方向网格数Z方向网格数X方向网格数X方向网格数强度ExV/M(0.015m)通过挡板两端都吸收(V)0.20.150.10.05第一个第二个观测面入射电压观测面透射电压第一个观测面反射电压-0.05J*EEE=EE时间步长010020030

15、04005006007008009001000参考面电压2.5*1011s)1终端匹配时S参数终端短路、入射口吸收,riEx被终端反射的波形z方向网格数(0.015m)频率f(0.1GHz)终端短路：由于挡板与短路板之间构成一种谐振腔结构，时域脉冲波在其间来回反射，入射口仍然是吸收边界条件，波完全消失的时间会比终端匹配情况要长。Nt足够大(Nt=20000)如下图所示:电场强度ExV/M电压幅度Vx10终端短路时第一参考面电压频率f(0.1GHz)反射系数(f)终端短路时反射系数频率f(0.1GHz)第五部分课程总结与讨论短短几周的电磁仿真综合实践很快就结束了，这几周紧张的学习给我留下了深刻的

16、印象。我们已经学习了电磁场与电磁波、微波技术基础等课程，对电磁理论有一定的认识和了解，但那只是完全理论性质的。通过这门课程的学习，我对电磁场等相关知识有了进一步的了解，有了一定感性的认识，感觉电磁场也不是那么和天书一样的难以理解。通过用计算机模拟将电磁现象展示在我们面前，使得我不但加深了对电磁理论的理解，同时也锻炼了我们计算机应用方面的能力，对科学研究有了一定的了解。我有如下几点体会：首先，这次实践使得我对计算机的应用有了进一步的提高。以往我对计算机的认识只是处于一个很简单的层次上，对其应用、原理也只是略知皮毛。但通过这段时间的学习，自己用计算机以及相关软件完成了这次实践活动。虽然在开始的时候

17、对软件以及原理理解不够，但是通过老师的多次讲解，同学之间的互相探讨，渐渐地自己掌握了一定的方法。看着自己编出的程序调试出令人满意的结果，还是相当喜悦的。其次，这次的实践活动让我加深了对专业知识的了解，对电磁理论也不象以前那样见到就躲了。当我用自己编的程序捕捉到了它们，对它们的产生、传播、反射、透射等现象有了更为直接的认识和了解。自己本来对电磁理论有一定的兴趣，苦于没有很好的方法帮助自己更好的理解和学习，或许通过这次实践活动，自己以后能够找到适合自己的学习电磁理论的方法。第三，通过对原理、程序的不断学习、操作，使得我认识到做任何事情都要认真，持之以恒。王老师科学严谨的作风给我留下了深刻的印象，不

18、论是做人还是做科研，王老师都是那么的一丝不苟。这一点让我也不敢有所懈怠。同时我也感觉到自己能力的欠缺，分析问题、解决问题的能力还有待提高。相信在以后我会通过自己的不断努力，提高自己的能力，学习和工作上更加细致认真。第四，这次实践活动，我感觉是我上大学以来收获最大的一门课程。王老师讲课不是传统的灌输式的教育，学生是真正意义上学习的主体。课前需要我们对理论方面的知识多加复习和学习，上课的时候需要我们积极的开动脑筋，如何把理论转化为实际可见，课后也需要我们不断的思考，如何才能使得我们做的软件更加趋于完善。在掌握计算电磁学有关理论知识的同时，更加重要的是我们学到了方法，这一点我觉得才是最为重要的。下面

19、说一下我对这门课程的一些建议吧。我觉得要消除一些学生对这门课程的紧张情绪（开始我就有这样的感觉），适当的增加课时数，毕竟我们目前涉及的真的只是皮毛而已。增加课时数，可以让大家对理论的理解加深，同时也会有更多的同学以后投入到电磁场领域方面的学习和研究。还有上机时间最好不要靠得太近了，毕竟学习有一个消化吸收的过程。这次的科研实践活动是令人愉快的。它让我体会到了快乐，当然也有抓耳挠腮的痛苦；它让我学习到了更多新的知识，新的方法；它让我看到了自己的不足之处，在以后的学习过程中，我一定会加以改正的。最后我还要特别感谢王老师，以及研究生师姐在这次实践活动中给予我的帮助。正是他们一丝不苟的工作态度和严谨的作

20、风深深的教育了我，以后我也会尽我最大的努力去做好每一件事情。第六部分附录一、参考文献1、计算电磁学(节选)»,王秉中，电子科技大学，2000年版2、电磁仿真综合实践，王秉中，电子科技大学3、MATLAB基础与应用简明教程，张平，北京航空航天大学出版社,2001年1月第一版二、相关程序1、一维终端匹配%定义变量t=0.5E-9;c=3E8;fmax=1/(2*t);lammin=c/fmax;epslon=8.854E-12;mu=4E-7*p;y0=120*pi;dz=lammin/20;dt=dz/(2*c);t0=3*t;l=6;d=0.18;nz=400;nt=1000;k1=

21、101;k2=301;cm=dt/(mu*dz);ce=dt/(epslon*dz);cmur=(dt*c-dz)/(dt*c+dz);ex(1:nz+1)=0;hy(1:nz)=0;vref1(1:nt)=0;vref2(1:nt)=0;讯序forn=1:ntex(1)=exp(-(n*dt-t0)A2/tA2);fori=1:nzhy(i)=hy(i)+cm*(ex(i)-ex(i+1);endexm=ex(nz);forj=2:nzex(j)=ex(j)+ce*(hy(j-1)-hy(j);end一维匹配终端吸收条件电场强度监控'g');pause(0.001);%观测面

22、电压图ex(nz+1)=exm+cmur*(ex(nz)-ex(nz+1);%vref1(n)=ex(k1)*d;vref2(n)=ex(k2)*d;%plot(1:nz+1,ex);pause(0.001);%x=1:nt;plot(x,vref1(x),'r',x,vref2(x),%gridend%一维终端匹配2、一维终端短路%定义变量y0=120*pi;dz=lammin/20;dt=dz/(2*c);t0=3*t;l=6;d=0.18;nz=400;nt=2000;k1=101;k2=301;cm=dt/(mu*dz);ce=dt/(epslon*dz);cmur=(

23、dt*c-dz)/(dt*c+dz);ex(1:nz+1)=0;hy(1:nz)=0;vref1(1:nt)=0;vref2(1:nt)=0;forn=1:ntfori=1:nzhy(i)=hy(i)+cm*(ex(i)-ex(i+1);endifn<=800ex(1)=exp(-(n*dt-t0)A2/tA2);forj=2:nzex(j)=ex(j)+ce*(hy(j-1)-hy(j);endelseexm=ex(2);forj=2:nzex(j)=ex(j)+ce*(hy(j-1)-hy(j);endex(nz+1)=0;%终端短路边界条件高斯波完全进入后，入口改为吸收边界条件&#

24、39;g');%两个观测面的电压波形实时监控Ex波形ex(1)=exm+cmur*(ex(2)-ex(1);%endvref1(n)=ex(k1)*d;vref2(n)=ex(k2)*d;x=1:nt;plot(x,vref1(x),'r',x,vref2(x),pause(0.01);%plot(1:nz+1,ex);pause(0.001);%gridend%一维终端短路3、二维终端匹配%定义变量t=0.5E-9;c=3E8;fmax=1/(2*t);lammin=c/fmax;epslon=8.854E-12;mu=4E-7*pi;y0=120*pi;dz=lam

25、min/20;dx=lammin/20;dt=dz/(2*c);t0=3*t;l=6;d=0.18;nz=400;nt=1000;k1=101;k2=301;cm=dt/(mu*dz);ce=dt/(epslon*dz);cf=dt/(epslon*dx);cmur=(dt*c-dz)/(dt*c+dz);nx=12;ez(1:nx+1,1:nz)=0;ex(1:nx,1:nz+1)=0;hy(1:nx,1:nz)=0;ver1in(1:nt)=0;ver1re(1:nt)=0;ver2tr(1:nt)=0;exb(1:nx,1)=0;exe(1:nx,1)=0;nf=10;df=1e8;jj

26、=sqrt(-1);gver1in(1:nf)=0;gver1re(1:nf)=0;gver2tr(1:nf)=0;s11(1:nf)=0;s21(1:nf)=0;forn=1:nt%时间循环fori=1:nx%二维FDTD过程fork=1:nzhy(i,k)=hy(i,k)+cm*(ez(i+1,k)-ez(i,k)*dz/dx+ex(i,k)-ex(i,k+1);endendifn<=300ex(1:12,1)=exp(-(n*dt-t0)A2/tA2);fori=1:nxfork=2:nzex(i,k)=ex(i,k)+ce*(hy(i,k-1)-hy(i,k);endendez(

27、1,1:nz)=0;ez(nx+1,1:nz)=0;ex(1:4,201)=0;ex(9:12,201)=0;%设置挡板fori=2:nxfork=1:nzez(i,k)=ez(i,k)+cf*(hy(i,k)-hy(i-1,k);endendez(1,1:nz)=0;ez(nx+1,1:nz)=0;ex(1:4,201)=0;ex(9:12,201)=0;else%高斯波完全进入，入射端改为吸收边界条件exb(1:nx,1)=ex(1:nx,2);exe(1:nx,1)=ex(1:nx,nz);fori=1:nxfork=2:nzex(i,k)=ex(i,k)+ce*(hy(i,k-1)-h

28、y(i,k);endendfori=2:nxfork=1:nzez(i,k)=ez(i,k)+cf*(hy(i,k)-hy(i-1,k);end入口边界条件endex(1:nx,1)=exb(1:nx,1)+cmur*(ex(1:nx,2)-ex(1:nx,1);%ex(1:nx,nz+1)=exe(1:nx,1)+cmur*(ex(1:nx,nz)-ex(1:nx,nz+1);%吸收边界条件ez(1,1:nz)=0;ez(nx+1,1:nz)=0;ex(1:4,201)=0;ex(9:12,201)=0;endfori=1:nx%分离出入射电压、反射电压、透射电压ver2tr(n)=ver2

29、tr(n)+ex(i,301)*d/12;endifn<=500fori=1:nxver1in(n)=ver1in(n)+ex(i,101)*d/12;endelsefori=1:nxver1re(n)=ver1re(n)+ex(i,101)*d/12;endend%x=1:nt%plot(x,ver1in(x),'r',x,ver1re(x),'b',x,ver2tr(x),'m');%各个电压波形surf(ex);%监控Ex波形pause(0.01)endforf=1:nfforn=1:nt%傅氏变换gver1in(f)=gver1in

30、(f)+dt*ver1in(n)*exp(-2*pi*f*df*n*dt*jj);gver1re(f)=gver1re(f)+dt*ver1re(n)*exp(-2*pi*f*df*n*dt*jj);gver2tr(f)=gver2tr(f)+dt*ver2tr(n)*exp(-2*pi*f*df*n*dt*jj);endendforf=1:nf%S参数s11(f)=gver1re(f)/gver1in(f);s21(f)=gver2tr(f)/gver1in(f);endfori=1:nfs1(i)=abs(s11(i);s2(i)=abs(s21(i);end%x=1:nf;%plot(x

31、,s1,'g',x,s2,'r');%S参数%grid%二维终端匹配4、二维终端短路%定义变量t=0.5E-9;c=3E8;fmax=1/(2*t);lammin=c/fmax;epslon=8.854E-12;mu=4E-7*pi;y0=120*pi;dz=lammin/20;dx=lammin/20;dt=dz/(2*c);t0=3*t;l=6;d=0.18;nz=400;nt=20000;k1=101;k2=301;cm=dt/(mu*dz);ce=dt/(epslon*dz);cf=dt/(epslon*dx);cmur=(dt*c-dz)/(dt*c+

32、dz);nx=12;ez(1:nx+1,1:nz)=0;ex(1:nx,1:nz+1)=0;hy(1:nx,1:nz)=0;vref1(1:nt)=0;vref2(1:nt)=0;ver1in(1:nt)=0;ver1re(1:nt)=0;ver2tr(1:nt)=0;exb(1:nx,1)=0;exe(1:nx,1)=0;nf=10;df=1e8;jj=sqrt(-1);gver1in(1:nf)=0;gver1re(1:nf)=0;gver2tr(1:nf)=0;s11(1:nf)=0;s21(1:nf)=0;forn=1:ntfori=1:nx%二维FDTD过程fork=1:nzhy(i,k)=hy(i,k)+cm*(ez(i+1,k)-ez(i,k)*dz/dx+ex(i,k)-ex(i,k+1);endendifn<=300ex(1:12,1)=exp(-(n*dt-t0)A2/tA2);fori=1:nxfork=2:nzex(i,k)=ex(i,k)+ce*(hy(i,k-1)-hy(i,k);endendez(1,1:nz)=0;ez(nx+1,1:nz)=0;ex(1:4,201)=0;ex(9:12,201)=0;