Sunday, January 18, 2009

Absorbing boundary conditions are necessary to keep outgoing and fields from being reflected back into the problem space. Normally, in calculating the field, we need to know the surrounding values; this is a fundamental assumption of the FDTD method. At the edge of the problem space we will not have the value to one side. However, we have an advantage because we know there are no sources outside the problem space. Therefore, the fields at the edge must be propagating outward. We will use these two facts to estimate the value at the end by using the value next to it.

Suppose we are looking for a boundary condition at the end where i=Imin. If a wave is going toward a boundary in free space, it is traveling at c0, the speed of light. So in one time step of the FDTD algorithm, it travels distance equals to derta/2. This equation basically explains that it takes two time steps for a wave front to cross one cell. So a common sense approach tells us that an acceptable boundary condition might be E(Imin,n)=E(Imin+1,n-2) and E(Imax,n)=E(Imax-1,n-2). This is the well know Mur absorbing boundary condition. It is relatively easy to implement as:


!* FD1D_002 1D FDTD simulation in free space *!
!* Absorbing Boundary Condition added
Module common_data
implicit none
!Calculation region
integer,save:: Imin,Imax
!constant
real,parameter:: c0=3.e8,Eps0=8.854e-12,Mu0=1.25663706e-6
!source parameter
real,save ::SourceT0,SourceT
!maximum number of time steps
integer,parameter:: timestop=500
!grid size, time step
real,save:: dx,dt
!variables for field
real,dimension(:),allocatable,save:: Ez,Hy
!variables for absorbing boundary
real,dimension(:),allocatable,save:: ABC
endmodule common _data

program main
use common _data
implicit none
integer ntime,i

call pre_processing
do Ntime=1,timestop
call Ezfld(Ntime)
call EzABC
call Hyfld
enddo
call savedata
end

subroutine pre_processing
use common _data
implicit none
Imin=-300;Imax=300
allocate(Ez(Imin:Imax))
allocate(Hy(Imin:Imax))
allocate(ABC(1:4))
Ez=0.;Hy=0.;ABC=0.
dx=1.e-2;
dt=dx/2./c0
SourceT0=120;SourceT=50
return
end

subroutine Ezfld(n)
use common _data
implicit none
integer i,n
do i=Imin+1,Imax-1
Ez(i)=Ez(i)+(Hy(i-1)-Hy(i))*dt/dx/Eps0
enddo
Ez(0)=Ez(0)+exp(-((n-1)*1.-SourceT0)**2/SourceT**2)
return
end

subroutine EzABC
use common _data
implicit none
Ez(Imin)=ABC(2)
ABC(2)=ABC(1)
ABC(1)=Ez(Imin+1)

Ez(Imax)=ABC(4)
ABC(4)=ABC(3)
ABC(3)=Ez(Imax-1)
return
end

subroutine Hyfld
use common _data
implicit none
integer i
do i=Imin,Imax-1
Hy(i)=Hy(i)-(Ez(i+1)-Ez(i))*dt/dx/Mu0
enddo
return
end
subroutine savedata
use common _data
implicit none
integer i
open(1,file='ex')
open(2,file='z_cell')
write(1,*) (Ez(i),i=Imin,Imax)
write(2,*) (i,i=Imin,Imax)
close(1)
close(2)
return
end


Fig.1 shows the results of a simulation using FD1D_002. A pulse that originates in the center propagates outward and is absorbed without reflecting anything back into the problem space.


No comments:

Post a Comment

Followers

Blog Archive