[讨论]球散射的外推程序，请大家帮忙看看

最近编了一段球的散射的程序，大家帮我看下我的这段频域三维远场外推的程序对不对 TL= Y
QA  
 ein0=cmplx(e0r,e0i)      原点处的入射电场 ydl jw  
 hin0=ein0/cexp(cmplx(0,0.5*del_t*2*pi*f))    入射电场滞后半个时间步 4kp im  
 do i=imin,imax D|I Ec?  
 do j=jmin,jmax ,zcQS-e2  
 do k=kmin,kmax k!jNOqbb  
 exf(i,j,k)=cmplx(exr(i,j,k),exi(i,j,k))/ein0           电场归一 UIJx*  
 eyf(i,j,k)=cmplx(eyr(i,j,k),eyi(i,j,k))/ein0 {hSGv
  
 ezf(i,j,k)=cmplx(ezr(i,j,k),ezi(i,j,k))/ein0 l2v_?j-)x  
 hxf(i,j,k)=cmplx(hxr(i,j,k),hxi(i,j,k))/hin0           磁场归一 Eu@
5L9A  
 hyf(i,j,k)=cmplx(hyr(i,j,k),hyi(i,j,k))/hin0 gLD`wfZR  
 hzf(i,j,k)=cmplx(hzr(i,j,k),hzi(i,j,k))/hin0 OQ/<-+<
w  
 end do 4}v@C|.p  
 end do !.-tW7  
 end do CDMfa&
;T  
       ~KjJ\b)R  
oqkVYlE  
c     外推边界的等效电磁流 lYf+V8{  
c     左平面 :2V^K&2L  
 do j=jomin,jomax-1 ]7sx;KFv  
 do k=komin,komax-1 $X`y%*<<v  
 jz(iomin,j,k)=-0.25*(hyf(iomin,j,k)+hyf(iomin,j+1,k) cl[!`Z  
     &     +hyf(iomin-1,j,k)+hyf(iomin-1,j+1,k)) :@ uIxa$[  
      jy(iomin,j,k)=0.25*(hzf(iomin,j,k)+hzf(iomin,j,k+1) R U"/2i  
     &     +hzf(iomin-1,j,k)+hzf(iomin-1,j,k+1)) n[xkSF^)  
 jmz(iomin,j,k)=0.5*(eyf(iomin,j,k)+eyf(iomin,j,k+1)) P=AS>N^yaL  
 jmy(iomin,j,k)=-0.5*(ezf(iomin,j,k)+eyf(iomin,j+1,k)) +%J\y^09kr  
 kr(iomin,j,k)=iomin*sin(thida)*cos(phi) 1N`vCt]w  
     &  +(j+0.5)*sin(thida)*sin(phi)+(k+0.5)*cos(phi) Ar9nBJ`  
 end do 6@;ha=[+  
 end do uE_c4Hp  
  LDT(]HJ  
c 右平面 KM?1/KZ/~  
      do j=jomin,jomax-1 %:6?Y%`*[  
 do k=komin,komax-1 /z."
l!u6  
 jz(iomax,j,k)=0.25*(hyf(iomax,j,k)+hyf(iomax,j+1,k) YZ#V#[j'^  
     &     +hyf(iomax-1,j,k)+hyf(iomax-1,j+1,k)) qcB){p+UQ  
      jy(iomax,j,k)=-0.25*(hzf(iomax,j,k)+hzf(iomax,j,k+1) S^RUw  
     &     +hzf(iomax-1,j,k)+hzf(iomax-1,j,k+1)) 3EFD%9n  
 jmz(iomax,j,k)=-0.5*(eyf(iomax,j,k)+eyf(iomax,j,k+1)) iMP  
 jmy(iomax,j,k)=0.5*(ezf(iomax,j,k)+eyf(iomax,j+1,k)) noVa=aU^  
 kr(iomax,j,k)=iomax*sin(thida)*cos(phi) `ue[q!Qq  
     &  +(j+0.5)*sin(thida)*sin(phi)+(k+0.5)*cos(phi) )yee2(S  
 end do S[M4ukYK  
 end do      Z%o7f6P0IX  
C     前平面 Lrz3  
      do i=iomin,iomax-1 '3tw<k!1{.  
 do k=komin,komax-1 Q}%tt=KD  
 jz(i,jomin,k)=0.25*(hxf(i,jomin,k)+hxf(i+1,jomin,k)  d<p2/aA  
     &     +hxf(i,jomin-1,k)+hxf(i+1,jomin-1,k)) O0l^*nZ46t  
 jx(i,jomin,k)=-0.25*(hzf(i,jomin,k)+hzf(i,jomin,k+1) /4S;QEv  
     &     +hzf(i,jomin-1,k)+hzf(i,jomin-1,k+1)) ^~ =9  
      jmz(i,jomin,k)=-0.5*(exf(i,jomin,k)+exf(i,jomin,k+1)) jdRq6U^  
      jmx(i,jomin,k)=0.5*(ezf(i,jomin,k)+ezf(i+1,jomin,k)) s}x>J8hK  
 kr(i,jomin,k)=(i+0.5)*sin(thida)*cos(phi) 8@K^|xeQ  
     & +jomin*sin(thida)*sin(phi)+(k+0.5)*cos(phi) mx
Tk+j=  
 end do 9 wa,k  
 end do %(m])  
C     后平面 q1Qje%9@t  
      do i=iomin,iomax-1 u/^|XOy  
 do k=komin,komax-1 }amU[U,  
 jz(i,jomax,k)=-0.25*(hxf(i,jomax,k)+hxf(i+1,jomax,k)  jrJR1npB  
     &     +hxf(i,jomax-1,k)+hxf(i+1,jomax-1,k)) Bl.u=I:Y4  
 jx(i,jomax,k)=0.25*(hzf(i,jomax,k)+hzf(i,jomax,k+1) <zfe}0  
     &     +hzf(i,jomax-1,k)+hzf(i,jomax-1,k+1)) oT+(W,G  
      jmz(i,jomax,k)=0.5*(exf(i,jomax,k)+exf(i,jomax,k+1)) ^%bBW6eZ  
      jmx(i,jomax,k)=-0.5*(ezf(i,jomax,k)+ezf(i+1,jomax,k)) 0Vv9BL{  
 kr(i,jomax,k)=(i+0.5)*sin(thida)*cos(phi) }Ag|gF!_  
     & +jomax*sin(thida)*sin(phi)+(k+0.5)*cos(phi) |g)FA_#|<  
 end do oVkq2  
 end do  Ng<1Sd|MV  
c     下平面 ~&G4)AM  
      do i=iomin,iomax-1 =WZ%H_oxi  
 do j=jomin,jomax-1 I@7/jUO  
 jx(i,j,komin)=0.25*(hyf(i,j,komin)+hyf(i,j+1,komin) <zB*'m  
     &     +hyf(i,j,komin-1)+hyf(i,j+1,komin-1)) j\
)H  
 jy(i,j,komin)=-0.25*(hxf(i,j,komin)+hxf(i+1,j,komin) ^tAO_~4  
     &     +hxf(i,j,komin-1)+hxf(i+1,j,komin-1)) z MdC  
      jmx(i,j,komin)=-0.5*(eyf(i,j,komin)+eyf(i+1,j,komin))   ,TC~~EWq  
 jmy(i,j,komin)=0.5*(exf(i,j,komin)+exf(i,j+1,komin)) rnhFqNT:  
 kr(i,j,komin)=(i+0.5)*sin(thida)*cos(phi) 6XUuGxQV/  
     & +(j+0.5)*sin(thida)*sin(phi)+komin*cos(phi) LVtu*k   
 end do +H'\3^C-  
 end do IhonnLLW  
c     上平面 a%
,fXp>  
      do i=iomin,iomax-1 9w^zY;Y  
 do j=jomin,jomax-1 3P=w =~e  
 jx(i,j,komax)=-0.25*(hyf(i,j,komax)+hyf(i,j+1,komax) =V$j6  
     &     +hyf(i,j,komax-1)+hyf(i,j+1,komax-1)) ksI>IW  
 jy(i,j,komax)=0.25*(hxf(i,j,komax)+hxf(i+1,j,komax) XF,<i1ZlM  
     &     +hxf(i,j,komax-1)+hxf(i+1,j,komax-1)) lglYJ,  
      jmx(i,j,komax)=0.5*(eyf(i,j,komax)+eyf(i+1,j,komax)) 2hOPzv&B  
 jmy(i,j,komax)=-0.5*(exf(i,j,komax)+exf(i,j+1,komax)) ~uaP$*B[  
 kr(i,j,komax)=(i+0.5)*sin(thida)*cos(phi) x!{   
     & +(j+0.5)*sin(thida)*sin(phi)+komax*cos(phi) -zfoRU v  
 end do @ RR\lZ  
 end do L@)&vn]  
      tB6k|cPC  
C 由于非外推边界面上的电磁流取为0，所以为了方便，取区域内所有电磁流的和 [B/0-(?
 
 do i=iomin,iomax ym%slg  
 do j=jomin,jomax us{nyil1  
 do k=komin,komax 1M5 -pZ[D  
 fx=fx+jx(i,j,k)*exp(cmplx(0,kr(i,j,k)*kk))        3个方向的电流矩，KK为波矢 PnWD}'0V  
      fy=fy+jy(i,j,k)*exp(cmplx(0,kr(i,j,k)*kk)) o"_=K%9  
      fz=fz+jz(i,j,k)*exp(cmplx(0,kr(i,j,k)*kk)) |qUrEGjiSS  
      fmx=fmx+jmx(i,j,k)*exp(cmplx(0,kr(i,j,k)*kk))       3个方向的磁流矩 w+UV"\!G)Q  
      fmy=fmy+jmy(i,j,k)*exp(cmplx(0,kr(i,j,k)*kk)) >]$aoA#  
      fmz=fmz+jmz(i,j,k)*exp(cmplx(0,kr(i,j,k)*kk)) 2"/yEg*=  
 end do (O5)wej   
 end do *Zkss   
 end do E&2tBrAq  
    5{l1A(b  
  m1=z0*(fx*cos(thida)*cos(phi)+fy*cos(thida)*sin(phi)          THIDA方向远场公式的一部分 2R@%Y/  
     &    -fz*sin(thida))+(-fmx*sin(phi)+fmy*cos(phi)) : }?{@#Z  
 t1=(abs(m1))**2      !Tfij(91  
                 _vrWj<wyf  
 m2=z0*(fx*sin(phi)-fy .. S ~|.&0"\  
9{eBg
dC  

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

还是想办法找个已知的例子验证吧。 }2xb&6g~o  
裸看code，有点痛苦。