FDTD(时域有限差分法)算法

% Program author: Susan C. Hagness

% Department of Electrical and Computer Engineering % University of Wisconsin-Madison

% 1415 Engineering Drive

% Madison, WI 53706-1691

% 608-265-5739

% hagness@

%

% Date of this version: February 2000

%

% This MATLAB M-file implements the finite-difference time-domain % solution of Maxwell's curl equations over a three-dimensional

% Cartesian space lattice comprised of uniform cubic grid cells.

%

% To illustrate the algorithm, an air-filled rectangular cavity

% resonator is modeled. The length, width, and height of the

% cavity are 10.0 cm (x-direction), 4.8 cm (y-direction), and

% 2.0 cm (z-direction), respectively.

% conditions:

% ex(i,j,k)=0 on the j=1, j=jb, k=1, and k=kb planes

% ey(i,j,k)=0 on the i=1, i=ib, k=1, and k=kb planes

% ez(i,j,k)=0 on the i=1, i=ib, j=1, and j=jb planes

% These PEC boundaries form the outer lossless walls of the cavity. %

% The cavity is excited by an additive current source oriented

% along the z-direction. The source waveform is a differentiated

% Gaussian pulse given by

% J(t)=-J0*(t-t0)*exp(-(t-t0)^2/tau^2),

% where tau=50 ps. The FWHM spectral bandwidth of this zero-dc- % content pulse is approximately 7 GHz. The grid resolution

% (dx = 2 mm) was chosen to provide at least 10 samples per

% wavelength up through 15 GHz.

%

% To execute this M-file, type "fdtd3D" at the MATLAB prompt.

% This M-file displays the FDTD-computed Ez fields at every other % time step, and records those frames in a movie matrix, M, which % is played at the end of the simulation using the "movie" command. %

%***********************************************************************

clear

%***********************************************************************

% Fundamental constants

%***********************************************************************

cc=2.99792458e8; %speed of light in free space

muz=4.0*pi*1.0e-7; %permeability of free space

epsz=1.0/(cc*cc*muz); %permittivity of free space

%*********************************************************************** % Grid parameters

%***********************************************************************

ie=50; %number of grid cells in x-direction

je=24; %number of grid cells in y-direction

ke=10; %number of grid cells in z-direction

ib=ie+1;

jb=je+1;

kb=ke+1;

is=26; %location of z-directed current source

js=13; %location of z-directed current source

kobs=5;

dx=0.002; %space increment of cubic lattice

dt=dx/(2.0*cc); %time step

nmax=500; %total number of time steps

%*********************************************************************** % Differentiated Gaussian pulse excitation

%***********************************************************************

rtau=50.0e-12;

tau=rtau/dt;

ndelay=3*tau;

srcconst=-dt*3.0e+11;

%*********************************************************************** % Material parameters

%***********************************************************************

eps=1.0;

sig=0.0;