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

我想在一维自由空间连接边界处沿+Z方向引入入射波，空间范围为（-1000：1000），连接边界设在Z=-100处，仿真结果入射波确实在Z=-100处产生了，但入射波沿正反方向传播，调了很久都调不出来，程序在下面，用fortran写的，请各位大侠帮帮忙，帮我看看问题出在哪？先谢谢了！ b;S6'7Jf9  
w Y8@1>ah  
a?5WKO  
X AQGG>  
Yo>`h2C4  
EkN_8(w  
module variant u?Fnlne4@  
   iI3v[S  
     parameter(za=-1000)                   !计算的空间范围 2]L=s3  
   parameter(zb=1000) an3~'g?  
   L +L9Y}  
      integer i,t                           !循环变量:空间循环,时间循环   iaY5JEV:CA  
   ,"~WkLI~\t  
   real Ex(za:zb)                        !定义场量 CPVzX%=  
   real Hy(za:zb) /_]ltXD  
      real Ein(-1000:1000) 1TKEm9j]u  
      real Hin(-1000:999) `_|aeoK_  
      real EBin(4) 6ayy[5tW  
      real dt,dz,pi,n,f                      !时间步长,空间步长,圆周率,总的时间步,信号频率 u3
U4UK  
   real c,e0,u0,ur                        !光速,介电常数,磁导率,相对磁导率 q_cC7p6t  
   real er(za:zb)                         !相对介电常数 xC<=~(  
   real t0 ,spread 3\;v5D:  
      real dzt                               !中间变量 hT?6sWa  
   
>BBl7  
   real Ex1(za:za+1)                      !处理吸收边界 AT"!{Y "H  
   real Ex2(zb-1:zb) 2zz,(RA  
ZTN(irK  
end module \ jECSV|  
 }VF#\q  
program main >xT^RYS  
   O/mR9[}  
   use variant 3)}(M
 
   implicit none GxxDY]!  
   9/#0
?(K8  
   call ini() Yl$X3wi  
538fK9[  
   do t=0,n FC6xFg^  
G3G/xC"  
      call source() +}P%HH]E/p  
      call caculate_H()                         Z+ixRch@-s  
      call conection_H ()   k!T|)\nc+  
      call caculate_er()                                       UcgG  
      call caculate_E() *P`v^&  
      call conection_E () djWcbC=g_  
      call mur() 2z4<N2!M  
      call putout2() ]fj-`==  
   end do VM7 !0  
`)T&~2n
 
   call putout() M&r2:Whk  
1n_;kaY  
end Q8p=!K  
F9-[%l  
jAXR`D  
subroutine ini()                            !初始化变量的值 K4C^m|e  
   X2qv^G,  
   use variant YydA6IK4  
   implicit none hR-K@fS%l'  
   aR _NyA  
   c=3.0e8 qP7G[%=v  
   dz=0.0012                                         UVrQV$g!  
   n=1200 %;7.9%  
   f=2.45e9 zxkM'8JC  
   pi=3.1415926 W;4Lkk$  
   e0=8.854e-12 9CB\n  
   u0=pi*4.0e-7 #<4--$Xo  
   er=1 ]fZ<`w8u}  
   t0=40 Jm4#V~w  
   spread=12 5y]io Jc9-  
   ur=1     -OrR $w|e  
                                               %Uy%kN_&  
   dzt=c/dz                                   !由COURANT稳定条件确定时间步   f-4.WW2FN  
  dt=1.0/sqrt(dzt*dzt) _K'7(d0z  
   Wbe0ZnM]  
   Ex=0                                          !初始场值为0 Jy]Id*u9  
   Hy=0 8!&nKy<Y  
   Ex1=0                                         !MUR边界储存变量初始为0 n@hl2M6.x9  
   Ex2=0 TjxA#D)   
joqWh!kv7U  
    Ein=0 Ry&q1j  
    Hin=0 )>\4ULR83  
    EBin=0 ,of]J|  
end ^-hErsK  
8'_ 0g[s  
subroutine caculate_er()                       !计算时变的相对介电常数 /t*YDWLg  
t0q_>T-kt  
   use variant {6h|6.S2  
   implicit none ND5$bq Nu?  
_-O cc=Z
 
   do i=20,zb iRNLKi  
      er(i)=1 Di>rO038  
   end do V>Fesm"aq  
eJHp6)2  
end B8H75sz  
P>(P2~$Y"  
subroutine conection_E ()                            !E连接边界 ?R|th Z  
   {pWBwf>R C  
   use variant nU2V]-qY  
   implicit none Z, T#,  
-p`L%xj\  
     ~M5:=zKQ  
      Ex(-100)=Ex(-100)+dt/dz/e0*Hin(-101) ; DDe.f"  
     R(G\wqHUT3  
end wO7t!35  
Mu{BUtkzG  
subroutine conection_H ()                             !H连接边界 dQD$K|aUp  
   D^+?|Y@N  
   use variant :Xfn@>;3ui  
   implicit none v>H=,.`0\  
RP! X8~8  
      Hy(-101)=Hy(-101)+dt/dz/e0*Ein(-100) $e }n  
     aE)by-'  
end }^Gd4[(,g  
?iv=53<c#  
8O7Yv<  
subroutine caculate_E()                               !电场递推 +6l]]*H  
   .IgRY\?Q  
   use variant px=]bALU  
   implicit none xr[Vp  
RF?DtNuq  
   do i=za+1,zb-1 gzf-)J  
      Ex(i)=Ex(i)-dt/dz/(e0*er(i))*(Hy(i)-Hy(i-1)) T  p<s1'"  
   end do wC`;f5->  
QB&BTT=!  
end RNWX.g)b  
 @pFj9[N  
AOb]qc  
subroutine caculate_H()                              !磁场递推 eh(<m8I  
     wfR&li{  
   use variant A_R!uRD8-  
   implicit none )K;]y-Us[  
8md*wEjk  
   do i=za+1,zb-1 Q9c)k{QZ  
      Hy(i)=Hy(i)-dt/dz/(u0*ur)*(Ex(i+1)-Ex(i)) !:3.D,  
   end do ;hHi@Z9  
end FJ54S
 
booRrTS  

F PR`tE  
subroutine source()                                 !入射波源 .WX,Nd3@  
   U]R~gy}#  
   use variant Q`ME@vz  
   implicit none I[d]!YI}F  
    do i=-1000,999 c?E{fD"Fc3  
        Hin(i)=Hin(i)-dt*(Ein(i+1)-Ein(i))/(dz*u0) P`(Mk6gE  
    end do `)Ky0&?  
    do i=-999,999 [=uI
b._Wv  
        Ein(i)=Ein(i)+dt*(Hin(i-1)-Hin(i))/(dz*e0)  &;c>O  
    end do $SY]fNJQ  
pK}=*y~$  
      Ein(-120)= 1*sin(2*pi*f*real(t)*dt)   y9L14  
     =-#G8L%Q  
    Ein(-1000)=EBin(4) fThgK;Qy'U  
      EBin(4)=EBin(3) V5p0h~PK  
    EBin(3)=Ein(-999) r2\c'9uH  
    Ein(1000)=EBin(2) [syj#  
    EBin(2)=EBin(1) _mJG5(|  
    EBin(1)=Ein(999) },@``&
e  
               HC*?DJ,  
end "& 25D  
,SUT~oETP  
subroutine mur() pHC/(6?  
   !L$x:/R9M  
   use variant !<<AzLVL  
   implicit none Ak@Dyi?p  
W0I4Vvh_"  
  Ex(za)=Ex1(za+1)+(c*dt-dz)/(c*dt+dz)*(Ex(za+1)-Ex1(za))  !入射波吸收边界条件(一阶MUR) H:QhrL+7_  
   Ex(zb)=Ex2(zb-1)+(c*dt-dz)/(c*dt+dz)*(Ex(zb-1)-Ex2(zb)) mJ5H=&Z  
     BPv>$ m+.  
   Ex1(za+1)=Ex(za+1) * F4UAQzYb  
   Ex1(za)=Ex(za) xZ.!d.rn  
   Ex2(zb-1)=Ex(zb-1)
6RoAl$}'  
   Ex2(zb)=Ex(zb) =qu(~]2(  
ru9zTZZD  
end /7"I#U^u
/  
     Vzg=@A#  
F<|t\KOW  
subroutine putout()                                     !输出 2C=Q8ayvX  
   N( Oyi  
   use variant (xpn`NA  
 &nbs .. O;uG?.\  
k N7Bd}  

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