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

我想在一维自由空间连接边界处沿+Z方向引入入射波，空间范围为（-1000：1000），连接边界设在Z=-100处，仿真结果入射波确实在Z=-100处产生了，但入射波沿正反方向传播，调了很久都调不出来，程序在下面，用fortran写的，请各位大侠帮帮忙，帮我看看问题出在哪？先谢谢了！ t4uxon
 
p4_uY7^6  
~h+3WuOv  
+-V?3fQ  
.?:~s8kB  
QV1%Zou  
module variant Us.jyg7_c  
   .JjuY'-Q  
     parameter(za=-1000)                   !计算的空间范围 aa]v7d  
   parameter(zb=1000) j pV  
   OPjh"Hv  
      integer i,t                           !循环变量:空间循环,时间循环   8G<.5!f7`N  
   
)}5rs  
   real Ex(za:zb)                        !定义场量 b=EZtk6>  
   real Hy(za:zb) `rQl{$9IC  
      real Ein(-1000:1000) }$U6lh/Ep  
      real Hin(-1000:999) Vd~k4  
      real EBin(4) 8=uljn/  
      real dt,dz,pi,n,f                      !时间步长,空间步长,圆周率,总的时间步,信号频率 Ne_>%P
|I_  
   real c,e0,u0,ur                        !光速,介电常数,磁导率,相对磁导率 Mq$=zsj  
   real er(za:zb)                         !相对介电常数 Vtri"G8 aB  
   real t0 ,spread SBt: `,  
      real dzt                               !中间变量 FFGTIT# {"  
   OD|&qsbL  
   real Ex1(za:za+1)                      !处理吸收边界 z"`?<A&u  
   real Ex2(zb-1:zb) fsu"Lc  
MLX.MUS  
end module u1<kdTxA N  
R?I3xb  
program main U5=J;[w}N  
   7]bqs"t  
   use variant %Kp}Wo6  
   implicit none p {3|W<  
   T94$}- 5/)  
   call ini() XzsK^E0R  
d0az#Yg!  
   do t=0,n b L~<~gA  
I*"]!z1  
      call source() x]T;W&s  
      call caculate_H()                         =QRZ(2Wq  
      call conection_H ()   P]GGnT
(!  
      call caculate_er()                                        VVY\W!  
      call caculate_E() { q<l]jn9  
      call conection_E () 0}N^l=jQ  
      call mur() ,;y^|X  
      call putout2() RJM(+5xQ|  
   end do oY:>pxSz<@  
cPSu!u}D  
   call putout() zFIKB9NUn  

hRu%> =7  
end m\=u/Zip  
+0DIN4Y(4  
HY:
n{=o  
subroutine ini()                            !初始化变量的值 50A_+f.7%  
   XMdYted  
   use variant uv!/DX#  
   implicit none \Qn8"I83AV  
   ^xrR3m*d  
   c=3.0e8 jL<.?HE  
   dz=0.0012                                         Ygb#U'|  
   n=1200 ~5b %~:  
   f=2.45e9 &$[{L)D  
   pi=3.1415926 OMU#Sx!6  
   e0=8.854e-12 uP
Pe"$  
   u0=pi*4.0e-7 $? Rod;  
   er=1 =%p{"<  
   t0=40 pV8,b  
   spread=12 }2S \-  
   ur=1     H Vy^^$  
                                               hV)I C9  
   dzt=c/dz                                   !由COURANT稳定条件确定时间步   'I2)-=ZL6  
  dt=1.0/sqrt(dzt*dzt) {JJ`|*H$_  
   W%o){+,  
   Ex=0                                          !初始场值为0 =k z;CS+  
   Hy=0 j06DP _9M  
   Ex1=0                                         !MUR边界储存变量初始为0 Oc>-jhx?  
   Ex2=0 4=Zlsp  
4w,}1uNEf  
    Ein=0 zoU.\]#C  
    Hin=0 !p$V7pFu6  
    EBin=0 "7DPsPs  
end 2_Qzc&"[ 4  
frt?*|:  
subroutine caculate_er()                       !计算时变的相对介电常数 _%]H}N Q  
]F~5l?4u#  
   use variant TH-^tw  
   implicit none EMo6$(  
jWz-7BO  
   do i=20,zb moGbBkO  
      er(i)=1 ~+F: QrXcI  
   end do 6dh PqL  
E:ytdaiT  
end FYJB.lAT  
~Y(M>u.+!  
subroutine conection_E ()                            !E连接边界 ~oI49Q&{  
   Dcf`+?3  
   use variant T>?~eYHXs  
   implicit none f W )  
v|xlI4  
     3^\y>
 
      Ex(-100)=Ex(-100)+dt/dz/e0*Hin(-101) BI s!  
     Q.Acmht#  
end  O>3'ylBQ  
&zR}jD>  
subroutine conection_H ()                             !H连接边界 -:QyWw/d  
   <}=
D?bXw  
   use variant 19EU[eb  
   implicit none \r2w@F{C  
U7W ct %  
      Hy(-101)=Hy(-101)+dt/dz/e0*Ein(-100) *=8)]_=f  
     W[jxfZD9v  
end 5^k#
fl2  

ocMf}"  
C"}x=cK  
subroutine caculate_E()                               !电场递推 U6Qeode  
   vl
D]!]V:h  
   use variant VXa]L4jJ9  
   implicit none :1iw_GhJf  
8Xr3q eh+  
   do i=za+1,zb-1 A28w/=e7  
      Ex(i)=Ex(i)-dt/dz/(e0*er(i))*(Hy(i)-Hy(i-1)) qhtc?A/0}  
   end do `C1LR,J  
$2.DZ  
end L%DL n  
Z6=!}a%  
"YBA$ef$  
subroutine caculate_H()                              !磁场递推 YANEdH`d  
     IO+z:D{  
   use variant U;31}'b  
   implicit none bMZ0%(q
 
gSC@uf  
   do i=za+1,zb-1 Pzqgg43Xf  
      Hy(i)=Hy(i)-dt/dz/(u0*ur)*(Ex(i+1)-Ex(i)) H`
P )  
   end do \@xnC$dd/  
end W)l&4#__(  
75F&s,4+  
mH?hzxa+  
subroutine source()                                 !入射波源 #O~Y[''C5X  
   }LRAe3N%8  
   use variant p#&6Ed*V  
   implicit none ,%C$~+xjM  
    do i=-1000,999 H\^^p!^)  
        Hin(i)=Hin(i)-dt*(Ein(i+1)-Ein(i))/(dz*u0) 4Vs;Y&t]  
    end do ?:ZH%R_`a  
    do i=-999,999 Zj JD@,j  
        Ein(i)=Ein(i)+dt*(Hin(i-1)-Hin(i))/(dz*e0)
s2E}+ #  
    end do \bg^E>-  
0MOAd!N  
      Ein(-120)= 1*sin(2*pi*f*real(t)*dt)   WM#!X!Vo  
     [Yv5Sw  
    Ein(-1000)=EBin(4) L*_xu _F  
      EBin(4)=EBin(3) ovl@[>OB  
    EBin(3)=Ein(-999) 8c~H![2u  
    Ein(1000)=EBin(2) S$#Awen"@  
    EBin(2)=EBin(1) oeIS&O.K  
    EBin(1)=Ein(999) Kq*^*vWC  
               9.R_=  
end MLDuo|?  
=FE,G*  
subroutine mur() <[q)2 5RL  
   V*?,r<(  
   use variant L/ZZe5I  
   implicit none ql4T@r3l}3  
OJ 2M_q)e  
  Ex(za)=Ex1(za+1)+(c*dt-dz)/(c*dt+dz)*(Ex(za+1)-Ex1(za))  !入射波吸收边界条件(一阶MUR) $qvNv[  
   Ex(zb)=Ex2(zb-1)+(c*dt-dz)/(c*dt+dz)*(Ex(zb-1)-Ex2(zb)) L]VK9qB  
     vD(;VeW[  
   Ex1(za+1)=Ex(za+1) RHxd6Gs"  
   Ex1(za)=Ex(za) ,ToED  
   Ex2(zb-1)=Ex(zb-1) vH=I#Ajar  
   Ex2(zb)=Ex(zb) ex['{
|a{  
y`N1I  
end a-W&/  
     ;8A_-
$  
"[%NXan  
subroutine putout()                                     !输出 j}|6k6t  
   ] =D+a&  
   use variant >A)he!I  
 &nbs .. V;z?m)ur  
emHi=[!i  

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