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

以下是我编写的二维TM波的FDTD 的FORTRAN程序 b\k]Jx  
是自由空间中的电磁波的传播模拟，没有其他的限制条件，只是想实现模拟电磁波在二维真空中的传播情况，但这个程序的结果却是不稳定的，就是电场随时间无限增大，我看了下这个程序是符合时间稳定条件的，现在怎么也找不到毛病出在哪里了，请高手帮忙看看问题出在哪里： diq}\'f  
o.+;]i}D  
没有加边界的，但我想只要时间步小点应该还不会发散，但不知道为什么还是不收敛 51#*8u+L  
;O"?6d0  
      program fdtd2D IGdiIhH~2  
real*8,dimension(:,:),allocatable::hx,hy,ez sVLvnX,  
real eps0,mu0,dt,dx,dy,c0,pi,f gq+SM i=  
integer pmin,pmax qpjY &3SI  
pi=3.14159265 t un}rdb  
c0=3.0e8  !光速 C"{^w
y{sL  
f=3.e12    !源频率 bjZJP\6  
dx=3.e-7  !300纳米 zC?'Qiuh*  
dy=dx          dGIdSQ~ _  
dt=dx/(2*c0)        l
& :EKh  
eps0=8.854e-12 "[bkdL<  
mu0=4.0e-7 'solCAy  
nn=100  !网格数 zA,vp^  
pmin=0 b+rxin".  
pmax=nn $@g]?*L:  
  open(1,file='D:\MSDEV\ex.txt') i]JTKL{\q  
allocate (hx(pmin:pmax,pmin:pmax)) /F*Y~>*% 1  
allocate (hy(pmin:pmax,pmin:pmax)) M&Uy42,MR  
allocate (ez(pmin:pmax,pmin:pmax)) <Vat@e  
do i=pmin,pmax 0Y
.z  
  do j=pmin,pmax o-,."|6  
    hx(i,j)=0. ?YhDjQs  
  end do o<`vh*U@,4  
end do p~pD`'%  
do i=pmin,pmax +>ju,;4WK  
  do j=pmin,pmax j{@O%fv=  
    hy(i,j)=0. ,
6)N.  
  end do z+"tAVB[i  
end do          q!Du J  
do i=pmin,pmax i! .]U@{k  
  do j=pmin,pmax 7#V7D6j1  
    ez(i,j)=0. Y"Cf84E  
  end do ,+=9Rp`md  
end do jUe@xis<T  
#EzBB*kP  
nt=20 !时间步数 j6 wFks  
_:x]'w%  
do n=1,nt =~D? K9o  
do j=pmin+1,pmax =o@;K~
-  
  if(j>=pmin+1) then pwFp<O"  
    do i=pmin+1,pmax ^_7|b[Bt  
          if (i.eq.50.and.j.eq.50) then 0YH5B5b  
                ez(i,j)=sin(2*pi*f*n*dt) !激励源        ! a)_3r]sv^  
          else ]?9*Vr:P^  
                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<I+]Hf9  
          end if 7CIN!vrC|1  
    end do uYJS=NGNA  
  end if BI'>\hX/V  
  if(j<=pmax-1) then UhxM85M;x  
      do i=pmin,pmax-1 "TZq")-  
              hx(i,j)=hx(i,j)-(dt/mu0)*(ez(i,j+1)-ez(i,j))/dy dkQ4D2W*\  
  &nbs .. v88v
r  
^Hd[+vAvR  

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

感谢“憨企鹅”，已经解决了 r- 8Awa  
^y+k6bE  
真空磁导系数少乘个Pi

有没有c的源程序啊