[RFEDA原创]新手学习FDTD，自己编写的程序

 本人是一新手，混迹于本论坛，从各位老大的帖子里面学到很多。 lmd0Q(I  
刚刚开始学习FDTD，这个是我自己编写的程序，希望可以对像我这样的新手有一定启发哈。 Mu:H'$"'H  
  h&Sl8$jVp  
程序题目要求见附件 [F>zM  
_}\&;  
clear  %设置参数 NR3IeTd  
mu=4*pi*10^(-7);%真空中的磁导率 g)"6|Z?D"  
eps=8.8542*10^(-12);% 真空介电常数 5<S1,u5  
T=20*10^(-12); q
aQ
 
t0=4*T; )/=J=xw2  
c=3*10^8; @.gCeMlOf  
f=1/T;%频率 dNyc|P`U  
longth=c/f;%波长 
\2LC
pN  
ds=longth/15;%网格分辨率 SnXLjJe  
dt=ds/(2*c);%时间步长 P_Zo}.{  
%吸收边界条件迭代方程中的系数 LRmO6>y  
a=(c*dt-ds)/(c*dt+ds); ^8\pJg_0  
aa=(c^2*eps*dt)/(2*(c*dt+ds)); ]~oM'?&!  
%角点迭代方程中的系数 InDR\=o  
bb=(c*dt-sqrt(2)*ds)/(c*dt+sqrt(2)*ds); Wqqo8Y~fq  
nt=500; vuK 5DG4  
nx=200; _%>.t  
ny=nx; xNONf4I:6J  
   fHgfI@{=j  
%---------------------------------------------------初始化--------------------------------------------------- ]A\n>Z!;  
        Ex=zeros(nx,ny); >Te{a*`"m:  
        Ey=zeros(nx,ny); !P;qc
 
        Hz=zeros(nx,ny); %:/;R_  
        Hz(100,100)=exp(-16); >y=%o~  
        A/.cNen  
        44]ae~@a  
        dz=zeros(nx,ny);% n-1时刻的Hz l%)XPb2$J  
        ex=zeros(nx,ny);% n-1/2时刻的Ex
Zaj<*?\  
        ey=zeros(nx,ny);% n-1/2时刻的Ey Tu"yoF  
    vD*KJ3(c  
        2<'gX>TW  
  5|/vc*m_0'  
%------------------------------------------------电场和磁场迭代---------------------------------------------------% ' ZB%McS  
for t=1:nt  %时间步数 ;J_d%  
   ~7k b4[  
  %磁场迭代 D +oo5  
    j@:LMR>  
          ,rN7X<s54  
          <&NR3^Eq  
           Ex(1:nx,1:ny-1)=ex(1:nx,1:ny-1)+(dt/(eps*ds))*(Hz(1:nx,2:ny)-Hz(1:nx,1:ny-1)); }%}eyLm(  
       qlzL<  
     k82'gJ;MC=  
           Ey(1:nx-1,1:ny)=ey(1:nx-1,1:ny)-(dt/(eps*ds))*(Hz(2:nx,1:ny)-Hz(1:nx-1,1:ny)); 9hq7:  
           gi`ZFq@  
           pulse=exp(-(((t*dt-t0)/T)^2));%激励源高斯脉冲 z4<h)hh"k6  
 %电场迭代 *7Mrng  
  F%Xq}LMd  
        (O&b:D/Y  
      for i=2:nx :r6 bw  
        for j=2:ny-1 bv] ZUF0  
            if (i==nx*1/2&j==ny*1/2) 6D+9f{~r  
            Hz(i,j)=dz(i,j)+(dt/(mu*ds))*(Ex(i,j)-Ex(i,j-1)-Ey(i,j)+Ey(i-1,j))+pulse; .=NK^  
            else d E
Xw=u  
            Hz(i,j)=dz(i,j)+(dt/(mu*ds))*(Ex(i,j)-Ex(i,j-1)-Ey(i,j)+Ey(i-1,j)); N0XGW_f  
            end LO2sP"9  
        end 9*x9sfCv9  
           &Y,Rm78  
      end  F8k1fmM]Y  
 %---------------------------------------------------设置Mur的吸收边界条件---------- isN"7y|r:X  
 %------------------------------% )Y:CV,`  
  leiza?[  
     for j=2:ny-1 WT0U)x( m5  
        Hz(1,j)=dz(2,j)+a*(Hz(2,j)-dz(1,j))+aa*(ex(1,j)-ex(1,j-1)+ex(2,j)-ex(2,j-1)); EEaFi8  
    end @tP,l$O&  
 %边界x=nx Zs4N0N{  
    for j=2:ny-1 ]l\J"*"aB  
       Hz(nx,j)=dz(nx-1,j)+a*(Hz(nx-1,j)-dz(nx,j))+aa*(ex(nx-1,j)-ex(nx-1,j-1)+ex(nx,j)-ex(nx,j-1)); abF_i#  
    end .d*vfE$  
  %边界y=1  4f1*?HX&  
    for i=2:nx-1 t EN%mK  
        Hz(i,1)=dz(i,2)+a*(Hz(i,2)-dz(i,1))-aa*(ey(i,1)+ey(i,2)-ey(i-1,2)-ey(i-1,1)); 9:m+mpL
=9  
    end p&q&Fr-  
  %边界y=ny oi!E v_h  
    for i=2:nx-1 [~;9Mi.XL  
        Hz(i,ny)=dz(i,ny-1)+a*(Hz(i,ny-1)-dz(i,ny))-aa*(ey(i,ny-1)-ey(i-1,ny-1)+ey(i,ny)-ey(i-1,ny)); 9K&b1O@Aj  
    end "44X'G8N  
  CR_A{(  
     [g/g(RL  
  %---------------------------------------------------------角点处理------------------------------------------------------------% qzH97<M}T  
    Hz(1,1)=dz(2,2)+bb*(Hz(2,2)-dz(1,1));%角点x=1,y=1 l>:?U  
   c1Hv^*Y  
    Hz(nx,1)=dz(nx-1,2)+bb*(Hz(nx-1,2)-dz(nx,1));%角点x=nx,y=1 t &ucqY  
     y
K{~  
    Hz(1,ny)=dz(2,ny-1)+bb*(Hz(2,ny-1)-dz(1,ny) .. T|=8jt,  
D]0#A|nF  

未注册仅能浏览部分内容，查看全部内容及附件请先 登录 或 注册

不错，支持

呵呵，多发帖子，多交流，进步会快很多的。

这个的确要支持下

同感……交流才能共同进步。

恩，贵在坚持啊，：）

不错，支持下

程序只是手段，结果才是目的。加油。

相互学习 共同进步啊  

赞扬~ !Z\Gv1  
 *$DD+]2  
学习中~