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

 本人是一新手，混迹于本论坛，从各位老大的帖子里面学到很多。 r'*$'QY-N  
刚刚开始学习FDTD，这个是我自己编写的程序，希望可以对像我这样的新手有一定启发哈。 N#ggT9>X  
  FLWVI4*  
程序题目要求见附件 !.6n=r8d  
,MXU]{  
clear  %设置参数 `D77CC]vU  
mu=4*pi*10^(-7);%真空中的磁导率 Y$ jX  
eps=8.8542*10^(-12);% 真空介电常数 C
<D$Y,[w  
T=20*10^(-12); `<nxXsLe  
t0=4*T; s/7Z.\  
c=3*10^8; h+W$\T)  
f=1/T;%频率 ov_l)vt  
longth=c/f;%波长 f'=u`*(b7  
ds=longth/15;%网格分辨率 r84^/+"T  
dt=ds/(2*c);%时间步长 3jAr"xc  
%吸收边界条件迭代方程中的系数 L?h?LZnq  
a=(c*dt-ds)/(c*dt+ds);
{mJ' Lb0;  
aa=(c^2*eps*dt)/(2*(c*dt+ds)); Y%?S:&GH  
%角点迭代方程中的系数 M(a%Qk?]/  
bb=(c*dt-sqrt(2)*ds)/(c*dt+sqrt(2)*ds); ~M}{rl.n=  
nt=500; Z5vdH5?!r  
nx=200; ~X3g_<b_8  
ny=nx; |SOLC  
   MZV$YD^S  
%---------------------------------------------------初始化--------------------------------------------------- ^m#tWb)f  
        Ex=zeros(nx,ny); =de'Yy:\-  
        Ey=zeros(nx,ny); 2Xe1qzvo  
        Hz=zeros(nx,ny); a$=~1@  
        Hz(100,100)=exp(-16); 5_;-Qw  
         T01Iu  
        eUy*
0  
        dz=zeros(nx,ny);% n-1时刻的Hz I0+6p8,  
        ex=zeros(nx,ny);% n-1/2时刻的Ex Nm"P8/-09  
        ey=zeros(nx,ny);% n-1/2时刻的Ey xF[%R{Mn'  
    !i
"zM}
 
        hoq2zDjD  
  :K6JrS  
%------------------------------------------------电场和磁场迭代---------------------------------------------------% zr|DC] 3  
for t=1:nt  %时间步数 76!LMNf  
   Rn?JMM]  
  %磁场迭代 8oJl ]  
    <MbhBIejr  
          R;!,(l  
          3&y-xZu]  
           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)); w^0hVrws=,  
       7+P;s,mi7  
     _
f^6F<!  
           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)); qnp}#BZ  
           Rf!v{\  
           pulse=exp(-(((t*dt-t0)/T)^2));%激励源高斯脉冲
Z<;W*6J  
 %电场迭代 H7Q$k4\l  
  /9pxEidVAS  
        ?h[HC"V/2  
      for i=2:nx (I~
  
        for j=2:ny-1 p}QDX*/sSu  
            if (i==nx*1/2&j==ny*1/2) p3A9<g  
            Hz(i,j)=dz(i,j)+(dt/(mu*ds))*(Ex(i,j)-Ex(i,j-1)-Ey(i,j)+Ey(i-1,j))+pulse; q\@Zf}  
            else m#'eDO
:  
            Hz(i,j)=dz(i,j)+(dt/(mu*ds))*(Ex(i,j)-Ex(i,j-1)-Ey(i,j)+Ey(i-1,j)); AIg4u(j  
            end ]Kt@F0U<o  
        end t1hQ0B  
           q}Q G<%VR  
      end  f8;?WSGyD2  
 %---------------------------------------------------设置Mur的吸收边界条件---------- X]=eC6M}:V  
 %------------------------------% MDF_Xr-hZ  
  d^=)n-!T  
     for j=2:ny-1 2K9X(th1  
        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)); Tdcc<T  
    end bny5e:= d  
 %边界x=nx "K(cDVQ  
    for j=2:ny-1 r]!#v{#.  
       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)); ^4fvV\ne_~  
    end E
<jajYj  
  %边界y=1  f]r*;YEc4  
    for i=2:nx-1 4:/V
|E\D  
        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));
94w)Yln  
    end !.'D"Me>  
  %边界y=ny ~]-n%J$q  
    for i=2:nx-1 !Vyf2xS"  
        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)); 'V?FeWp  
    end ?>B?*IK!  
  t"4* ]S  
    Mc,|C)  
  %---------------------------------------------------------角点处理------------------------------------------------------------% 1b3 a(^^E  
    Hz(1,1)=dz(2,2)+bb*(Hz(2,2)-dz(1,1));%角点x=1,y=1 ?z.?(xZ 6  
   .Exvuo`F  
    Hz(nx,1)=dz(nx-1,2)+bb*(Hz(nx-1,2)-dz(nx,1));%角点x=nx,y=1 <YW)8J  
     @%,~5{Ir  
    Hz(1,ny)=dz(2,ny-1)+bb*(Hz(2,ny-1)-dz(1,ny) .. Fl_}Auj{&(  
qTJ0}F  

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

不错，支持

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

这个的确要支持下

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

恩，贵在坚持啊，：）

不错，支持下

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

相互学习 共同进步啊  

赞扬~ Vl?R?K=`~J  
d;(L@9HHD  
学习中~