[资料共享]看看我编的二维FDTD程序

以下是我编写的二维TM波的FDTD 的FORTRAN程序 LBTf}T\  
是自由空间中的电磁波的传播模拟，没有其他的限制条件，只是想实现模拟电磁波在二维真空中的传播情况，但这个程序的结果却是不稳定的，就是电场随时间无限增大，我看了下这个程序是符合时间稳定条件的，现在怎么也找不到毛病出在哪里了，请高手帮忙看看问题出在哪里： |(]XZ!{  
)Zox;}WK+  
没有加边界的，但我想只要时间步小点应该还不会发散，但不知道为什么还是不收敛 k/>k&^?  
W 5-=,t  
      program fdtd2D wtK+\Qnb  
real*8,dimension(:,:),allocatable::hx,hy,ez ~OXPn9qPp  
real eps0,mu0,dt,dx,dy,c0,pi,f Zv-#v  
integer pmin,pmax Vf0m7BJc3  
pi=3.14159265 @PutUYz  
c0=3.0e8  !光速 *5u0`k^j  
f=3.e12    !源频率 {>9vm!<[*\  
dx=3.e-7  !300纳米 r+>gIX+Fl  
dy=dx          n}EH{k9#  
dt=dx/(2*c0)        T)MKhK9\Ab  
eps0=8.854e-12 cJMp`DQzc  
mu0=4.0e-7 ?
g}G#j  
nn=100  !网格数 `KB;3L  
pmin=0 S3$&}I <  
pmax=nn %AaZc=a[c  
  open(1,file='D:\MSDEV\ex.txt') b-,4< H8m  
allocate (hx(pmin:pmax,pmin:pmax)) K9]L>Wj  
allocate (hy(pmin:pmax,pmin:pmax)) UZz/v#y~  
allocate (ez(pmin:pmax,pmin:pmax)) ;O+= 6>W  
do i=pmin,pmax 0 *2^joUv  
  do j=pmin,pmax F P3{Rp  
    hx(i,j)=0. !ap}+_IA7^  
  end do G8}k9?26(  
end do <v%
Q|r  
do i=pmin,pmax G .NGS%v  
  do j=pmin,pmax dSkMA  
    hy(i,j)=0. HQX.oW  
  end do KSz;D+L\  
end do          3s;^p,9 Y  
do i=pmin,pmax \8g'v@$wG  
  do j=pmin,pmax ew?4;  
    ez(i,j)=0. B 1jeIk,  
  end do #A\@)wJ  
end do $Y,y~4I  
evjj~xkt
e  
nt=20 !时间步数 .~u[rc|<  
DHQS7%)f`  
do n=1,nt fN&@y$  
do j=pmin+1,pmax JVydTvc  
  if(j>=pmin+1) then )Vd^#p  
    do i=pmin+1,pmax a`I \19p]  
          if (i.eq.50.and.j.eq.50) then Ym%XCl  
                ez(i,j)=sin(2*pi*f*n*dt) !激励源        ! TXdo,DPv7  
          else 4K5  
                ez(i,j)=ez(i,j)+(dt/eps0)*((hy(i,j)-hy(i-1,j))/dx-(hx(i,j)-hx(i,j-1))/dy) {>g{+Eq  
          end if rny(8z%Ck-  
    end do %tB7 &%ut  
  end if |yS  %  
  if(j<=pmax-1) then > 9.%hSy  
      do i=pmin,pmax-1 KrdEB0qh  
              hx(i,j)=hx(i,j)-(dt/mu0)*(ez(i,j+1)-ez(i,j))/dy B`|f"+.  
  &nbs .. A*G ~#v^  
7>"dc+Fg  

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

感谢“憨企鹅”，已经解决了 (n`] sbx  
9 ,:#Q<UM  
真空磁导系数少乘个Pi

有没有c的源程序啊