1379阅读
0回复

[求助]求助各位大侠：有关一维自由空间中连接边界处引入入射波问题

UID ：21885

注册：2008-11-23
登录：2014-10-29
发帖：33
等级：仿真新人

0楼发表于: 2009-12-16 22:06:50

我想在一维自由空间连接边界处沿+Z方向引入入射波，空间范围为（-1000：1000），连接边界设在Z=-100处，仿真结果入射波确实在Z=-100处产生了，但入射波沿正反方向传播，调了很久都调不出来，程序在下面，用fortran写的，请各位大侠帮帮忙，帮我看看问题出在哪？先谢谢了！

module variant

     parameter(za=-1000)                   !计算的空间范围
   parameter(zb=1000)

      integer i,t                           !循环变量:空间循环,时间循环

   real Ex(za:zb)                        !定义场量
   real Hy(za:zb)
      real Ein(-1000:1000)
      real Hin(-1000:999)
      real EBin(4)
      real dt,dz,pi,n,f                      !时间步长,空间步长,圆周率,总的时间步,信号频率
   real c,e0,u0,ur                        !光速,介电常数,磁导率,相对磁导率
   real er(za:zb)                         !相对介电常数
   real t0 ,spread
      real dzt                               !中间变量

   real Ex1(za:za+1)                      !处理吸收边界
   real Ex2(zb-1:zb)

end module

program main

   use variant
   implicit none

   call ini()

   do t=0,n

      call source()
      call caculate_H()
      call conection_H ()
      call caculate_er()
      call caculate_E()
      call conection_E ()
      call mur()
      call putout2()
   end do

   call putout()

end

subroutine ini()                            !初始化变量的值

   use variant
   implicit none

   c=3.0e8
   dz=0.0012
   n=1200
   f=2.45e9
   pi=3.1415926
   e0=8.854e-12
   u0=pi*4.0e-7
   er=1
   t0=40
   spread=12
   ur=1

   dzt=c/dz                                   !由COURANT稳定条件确定时间步
  dt=1.0/sqrt(dzt*dzt)

   Ex=0                                          !初始场值为0
   Hy=0
   Ex1=0                                         !MUR边界储存变量初始为0
   Ex2=0

    Ein=0
    Hin=0
    EBin=0
end

subroutine caculate_er()                       !计算时变的相对介电常数

   use variant
   implicit none

   do i=20,zb
      er(i)=1
   end do

end

subroutine conection_E ()                            !E连接边界

   use variant
   implicit none


      Ex(-100)=Ex(-100)+dt/dz/e0*Hin(-101)

end

subroutine conection_H ()                             !H连接边界

   use variant
   implicit none

      Hy(-101)=Hy(-101)+dt/dz/e0*Ein(-100)

end

subroutine caculate_E()                               !电场递推

   use variant
   implicit none

   do i=za+1,zb-1
      Ex(i)=Ex(i)-dt/dz/(e0*er(i))*(Hy(i)-Hy(i-1))
   end do

end

subroutine caculate_H()                              !磁场递推

   use variant
   implicit none

   do i=za+1,zb-1
      Hy(i)=Hy(i)-dt/dz/(u0*ur)*(Ex(i+1)-Ex(i))
   end do
end

subroutine source()                                 !入射波源

   use variant
   implicit none
    do i=-1000,999
        Hin(i)=Hin(i)-dt*(Ein(i+1)-Ein(i))/(dz*u0)
    end do
    do i=-999,999
        Ein(i)=Ein(i)+dt*(Hin(i-1)-Hin(i))/(dz*e0)
    end do

      Ein(-120)= 1*sin(2*pi*f*real(t)*dt)

    Ein(-1000)=EBin(4)
      EBin(4)=EBin(3)
    EBin(3)=Ein(-999)
    Ein(1000)=EBin(2)
    EBin(2)=EBin(1)
    EBin(1)=Ein(999)

end

subroutine mur()

   use variant
   implicit none

  Ex(za)=Ex1(za+1)+(c*dt-dz)/(c*dt+dz)*(Ex(za+1)-Ex1(za))  !入射波吸收边界条件(一阶MUR)
   Ex(zb)=Ex2(zb-1)+(c*dt-dz)/(c*dt+dz)*(Ex(zb-1)-Ex2(zb))

   Ex1(za+1)=Ex(za+1)
   Ex1(za)=Ex(za)
   Ex2(zb-1)=Ex(zb-1)
   Ex2(zb)=Ex(zb)

end


subroutine putout()                                     !输出

   use variant
&nbs ..

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

发帖回复