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

 本人是一新手，混迹于本论坛，从各位老大的帖子里面学到很多。 :4D#hOI  
刚刚开始学习FDTD，这个是我自己编写的程序，希望可以对像我这样的新手有一定启发哈。 'k9hzk(*  
  y'k4>,`9e  
程序题目要求见附件 ?JD\pYg[/  
\..(!>,%F  
clear  %设置参数 IBzHXa>75  
mu=4*pi*10^(-7);%真空中的磁导率 |KFWW  
eps=8.8542*10^(-12);% 真空介电常数 Px4zI9;cB  
T=20*10^(-12); r4c3t,L*$I  
t0=4*T; h;KI2k_^  
c=3*10^8; r4xq%hy  
f=1/T;%频率 rS 4'@a  
longth=c/f;%波长 ]&ptld;  
ds=longth/15;%网格分辨率 6]}Xi:I
 
dt=ds/(2*c);%时间步长 g/q$;cB  
%吸收边界条件迭代方程中的系数 ku'%+svD  
a=(c*dt-ds)/(c*dt+ds); []Z| *+=Q  
aa=(c^2*eps*dt)/(2*(c*dt+ds)); (;T;?v`-  
%角点迭代方程中的系数 `<>8tZS9"  
bb=(c*dt-sqrt(2)*ds)/(c*dt+sqrt(2)*ds); +Hb6j02#  
nt=500; ~QsQ7SAs  
nx=200; "A:wWb<
m  
ny=nx; GSW{h[Op  
   o+
O\VNW  
%---------------------------------------------------初始化--------------------------------------------------- %KY&E>^  
        Ex=zeros(nx,ny); 9&W\BQ  
        Ey=zeros(nx,ny); 'eo2a&S2D  
        Hz=zeros(nx,ny); ]~(Ipz2NP  
        Hz(100,100)=exp(-16); ^aZ Wu|p  
        7q|(ZZa  
        )vS##-[_  
        dz=zeros(nx,ny);% n-1时刻的Hz -c=IO(B/  
        ex=zeros(nx,ny);% n-1/2时刻的Ex DrltxI)  
        ey=zeros(nx,ny);% n-1/2时刻的Ey 1{"fmV  
    kkrQ;i)
Z  
        ;<q@>p[  
  zoOm[X=?3  
%------------------------------------------------电场和磁场迭代---------------------------------------------------% 5Tpn`2F  
for t=1:nt  %时间步数 ^%d+nKx9nL  
   ){>;eky  
  %磁场迭代 ]'_z(s}  
    !T0IMI  
          kQIWDN  
          V_}`
2.Pg  
           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)); oUw-l_M]  
       SQ5*?u\  
     X2i<2N*@  
           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)); U)p2PTfB  
           (p
14{  
           pulse=exp(-(((t*dt-t0)/T)^2));%激励源高斯脉冲 %e%nsj6  
 %电场迭代 o D* '  
  BJj~fNm1Zr  
        Y=3Y~  
      for i=2:nx 8d)F#  
        for j=2:ny-1 l.BSZhO$  
            if (i==nx*1/2&j==ny*1/2) ooZ7HTP|  
            Hz(i,j)=dz(i,j)+(dt/(mu*ds))*(Ex(i,j)-Ex(i,j-1)-Ey(i,j)+Ey(i-1,j))+pulse; GN1cnM>`  
            else Ky nZzR  
            Hz(i,j)=dz(i,j)+(dt/(mu*ds))*(Ex(i,j)-Ex(i,j-1)-Ey(i,j)+Ey(i-1,j)); 5Ll[vBW  
            end &7DE$ S  
        end $;;?'!%.  
           RL&3 P@r  
      end  0KDDAkR5R  
 %---------------------------------------------------设置Mur的吸收边界条件---------- bY>o%L
L-  
 %------------------------------% 6PMu;#  
  pb{P[-f  
     for j=2:ny-1 N{9<Tf*  
        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)); OC>" +  
    end ]owH [wvX  
 %边界x=nx C6P6hJm  
    for j=2:ny-1 sG VC+!E  
       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)); e8lF$[i  
    end . t3@86xTJ  
  %边界y=1  mtOrb9`m  
    for i=2:nx-1 {9l4 pT3  
        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)); 84knoC  
    end k2@IJ~  
  %边界y=ny z{n=G  
    for i=2:nx-1 yQx>h6  
        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)); vGAPQg6*  
    end W.67, 0m$  
  sJ?kp^!g  
    1Xs!ew)>  
  %---------------------------------------------------------角点处理------------------------------------------------------------% lT:<ZQyjT  
    Hz(1,1)=dz(2,2)+bb*(Hz(2,2)-dz(1,1));%角点x=1,y=1 Z9sg6M@s  
   {[9^@k  
    Hz(nx,1)=dz(nx-1,2)+bb*(Hz(nx-1,2)-dz(nx,1));%角点x=nx,y=1 Sc0ZT/Lm  
     fzKKK+  
    Hz(1,ny)=dz(2,ny-1)+bb*(Hz(2,ny-1)-dz(1,ny) .. A!f0AEA,  
Rxli;blzi  

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

不错，支持

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

这个的确要支持下

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

恩，贵在坚持啊，：）

不错，支持下

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

相互学习 共同进步啊  

赞扬~ L*v93;|s  
v/x*]c!"`  
学习中~