1、二维 TE 波 FDTD 法 PML 吸收边界 摘要:由于计算空间有限,一段时间后波面会被反射回来,这是我们不愿看到的。为了避免反射波的存在,需要在计算边界设置适当的吸收边界。目前最为有效的吸收边界是 PML 层。关键词:二维; TE 波;PML 层1 引言时域有限差分(F D T D )方法已逐渐成为解决电磁散射和电磁波传输等问题的有力工具。研究 F D T D 方法的核心问题是寻求一种理想的吸收边界, 使截断面反射最小。1 9 9 4 年 J.P.Berenger 提出了“ 完全匹配层(P M L ) ”这种新边界完全匹配层技术在计算电磁学领域特别是在 FEM 中得到广泛的应用是在 Sac
2、ks995 年提出的单轴各向异性介质 PLM 以后的事,此后的发展证明,作为一种新型的边界截断策略,不管是与 FDTD 的结合,还是用在 FEM 中,PLM 技术都获得了巨大成功。2 计算模型有一个长和宽都为 ke=80 的矩形介质,在它的周围围上厚度为 kp=8 的 PLM 层吸收边界(如图一所示) ,请编写程序验证电磁波在经过 PLM 层时完全被吸收(即电磁波没有反射)。图一模型图图二左为程序运行结果 Ey 图右为 Hzy 图图三为程序运行结果 hz 图3 结论:本文的算例证明了 PLM 作为一种边界处理的有效性,PLM 技术也在有关键参数的选取优化问题,这需要研究大量的算例来推出其规律,
3、因而还有很多的工作有待研究。参考文献:1z S Sacks,D M Kingsland,R Lee and Jin-Fa LeeAperfectly matched anisotropic absorber for use as anabsorbing boundary conditionJIEEE TransAntennasPropagat ,1995,43(12):146014632A Mitchell,T Aberle,D M Kokotoff and M WAustimAn anisotropic PML for use with biaxial mediaDIEEE TransMic
4、rowave Theory Tech ,1999,47(3):3743773班永灵填充各向异性介质的背腔式微带天线的矢量有限元法分析D北大学硕士论文,2003 年 5 月4班永灵(1978 一),男,河南人,电子科技大学博士研究生,主要研究方向有限元法(特剐是高阶矢量有限元)、快速扫频技术以及并行计算。5周乐柱(1944 一),男,贵州人,北京大学教授、博士生导师,CIE 会士,IEEE 高级会员,主要研究方向:计算电磁学及其应用(散射、天线、微波器件),通信中的电磁场问题等。6聂在平(1946 一),男,陕西人,电子科技大学教授,博士生导师,现任电子科技大学副校长,主要研究方向:计算电磁学、
5、电磁散射与逆散射、非均匀介质中的场与渡、移动通信中智能天线技术等。附录:二维 TE 波 FDTD 法 PML 吸收边界程序clc,clearke=80;kp=8;ex=ones(ke+2*kp+1,ke+2*kp+1);e=ones(ke+2*kp+1,ke+2*kp+1);ey=ones(ke+2*kp+1,ke+2*kp+1);hz=ones(ke+2*kp+1,ke+2*kp+1);hzx=ones(ke+2*kp+1,ke+2*kp+1);hzy=ones(ke+2*kp+1,ke+2*kp+1);kc=(ke+2*kp)/2;e0=8.85*10-12;mu0=4*3.14*10-7
6、;pi=3.14;c0=3e8;f=1e9;lambda=c0/f;dx=lambda/10;dy=lambda/10;dt=(dx/2)/c0;sinxm=0.333;sinmxm=sigxm*(mu0/e0);nn=3;T=0;for n=1:1000;T=T+1;hz=hzx+hzy;e=ex+ey;hz(kc,kc)=4*sin(2*pi*f*T*dt);for i=1:(2*kp+ke+1);for j=(kp+1):(ke+kp+1);ex(i,j)=ex(i,j)+dt/(e0*dx)*(hz(i,j)-hz(i,j-1);endendfor i=1:(ke+2*kp+1);fo
7、r j=2:kp;sinx=sinxm*(kp+2-j)/kp)nn*2*e0/dt;ex(i,j)=ex(i,j)*exp(-sigx*dt/e0)+(1-exp(-sigx*dt/e0)/(dy*sigx)*sqrt(e0/mu0)*(hz(i,j)-hz(i,j-1);endendfor i=1:(ke+2*kp+1);for j=(ke+kp+2):(ke+2*kp+1);sinx=sinxm*(j-ke-kp)/kp)nn*2*e0/dt;ex(i,j)=ex(i,j)*exp(-sigx*dt/e0)+(1-exp(-sigx*dt/e0)/(dy*sigx)*sqrt(e0/mu
8、0)*(hz(i,j)-hz(i,j-1);endendfor i=(kp+1):(ke+kp+1);for j=1:(2*kp+ke+1);ey(i,j)=ey(i,j)-dt/(e0*dx)*(hz(i,j)-hz(i-1,j);endendfor i=2:kp;for j=1:(ke+2*kp+1);sinx=sinxm*(kp+1-i)/kp)nn*2*e0/dt;ey(i,j)=ey(i,j)*exp(-sigx*dt/e0)-(1-exp(-sigx*dt/e0)/(dx*sigx)*sqrt(e0/mu0)*(hz(i,j)-hz(i-1,j);endendfor i=(ke+k
9、p+2):(ke+2*kp+1);for j=1:(ke+2*kp+1);sinx=sinxm*(i-ke-kp-1)/kp)nn*2*e0/dt;ey(i,j)=ey(i,j)*exp(-sinx*dt/e0)-(1-exp(-sigx*dt/e0)/(dx*sigx)*sqrt(e0/mu0)*(hz(i,j)-hz(i-1,j);endend%磁场循环开始for i=1:(ke+2*kp+1);for j=1:kp;sin(mx)=sinmxm*(kp-j+1)/kp)nn*2*mu0/dt;hzy(i,j)=hzy(i,j)*exp(-sigmx*dt/mu0)+(1-exp(-sig
10、mx*dt/mu0)/(dx*sigmx)*sqrt(mu0/e0)*(ex(i,j+1)-ex(i,j);endendfor i=1:(2*kp+ke+1);for j=(kp+1):(kp+ke);hzy(i,j)=hzy(i,j)+dt/(mu0*dy)*(ex(i,j+1)-ex(i,j);endendfor i=1:(ke+2*kp+1);for j=(kp+ke+1):(ke+2*kp);sigmx=sigmxm*(j-kp-ke)/kp)nn*2*mu0/dt;hzy(i,j)=hzy(i,j)*exp(-sigmx*dt/mu0)+(1-exp(-sigmx*dt/mu0)/(
11、dx*sigmx)*sqrt(mu0/e0)*(ex(i,j+1)-ex(i,j);endendfor i=1:kp;for j=1:(ke+2*kp+1);sigmx=sigmxm*(kp-i+1)/kp)nn*2*mu0/dt;hzx(i,j)=hzx(i,j)*exp(-sigmx*dt/mu0)-(1-exp(-sigmx*dt/mu0)/(dx*sigmx)*sqrt(mu0/e0)*(ey(i+1,j)-ey(i,j);endendfor i=(kp+1):(kp+ke);for j=1:(2*kp+ke+1);hzx(i,j)=hzx(i,j)+dt/(mu0*dy)*(ey(i
12、+1,j)-ey(i,j);endendfor i=(ke+kp+1):(ke+2*kp);for j=1:(ke+2*kp+1);sinmx=sinmxm*(i-kp-ke)/kp)nn*2*mu0/dt;hzx(i,j)=hzx(i,j)*exp(-sigmx*dt/mu0)-(1-exp(-sigmx*dt/mu0)/(dx*sigmx)*sqrt(mu0/e0)*(ey(i+1,j)-ey(i,j);endendey(1,:)=0;ex(:,1)=0;hzx(2*kp+ke+1,:)=0;hzy(:,2*kp+ke+1)=0;plot(hz(:,kc)axis(46,ke+2*kp,45,50)figure(1)subplot(1,2,1)mesh(ey)drawnowsubplot(1,2,2)mesh(hzy)drawnowfigure(2)mesh(hz)drawnowendhz%pcolor(hz)figure(3)contour3(hz,5)axis(kp,ke+kp,kp,kp+ke)shading interp%print-djpeg D2pml 相位图%print-djepg D2pml 等高线图