parameters.h

该代码可以用来计算二维TM波的FDTD程序。 外加理想匹配层（PML）。 有一定的广泛性。

#include<stdio.h>
#include<math.h>
#include<float.h>
#include<stdlib.h>
#include<conio.h>
#include<iomanip.h>
#include<afx.h>
#include<iostream.h>
#include<fstream.h>
#include<memory.h>

const double	pii=3.1415926535898;
const double	Epsilon0=1e-9/(36.*pii);	        //Permit in free space
const double	Mu0=1e-7*4.*pii;			//Permeab in free space
const double	LightSpeed=3.0e+8;			//Light Speed in free space(m/s)
const double	Eta0=120*pii;				//wave impendance in free space(Omn)
const double	Epsilon_r=1.0;				//Epsilon/Epsilon0
const double	Mu_r=1.0;				//Mu/Mu0
const double	Epsilon=Epsilon_r*Epsilon0;		//air space
const double	Mu=Mu_r*Mu0;

const double	Frequency=3.0e9;		        //frequency(3GHz)
const double	Wavelength=LightSpeed/Frequency;	//wavelength(0.1m)

//for fdtd: 2 dimension cartesian coordinates
const double	dx=Wavelength/140.0;//2.5e-4;
const double	dy=Wavelength/140.0;//2.5e-4;//
const double	dt=0.5*__min(dx,dy)/LightSpeed;//dx,dy=0.000714m
const double	tao=6e-11/dt;//50.4		

const int		Target_X=18;			//the X/2 thickness of scattering target　总场区x方向一半的长度
const int		Target_Y=18;			//the Y/2 thickness of scattering target
const int		rPlasma=0;//70;			//the plasma thickness 等离子区
const int		PML_Thickness=5;		//the thickness of PML
const int		Interval_PML_scatter=5;	//the interval between PML(吸收边界) and total-scattered field boundary(输出边界) 
const int		Interval_scatter_target=5;	//the interval between total-scattered field boundary(输出边界) and target(连接边界)

const int		MX=2*(Target_X+PML_Thickness+Interval_scatter_target+Interval_PML_scatter+rPlasma); //2*(18+5+5+5)=66
const int		MY=2*(Target_Y+PML_Thickness+Interval_scatter_target+Interval_PML_scatter+rPlasma);	//=66(in deltas)		
const int		dia=int(sqrt(MX*MX+MY*MY))+2*PML_Thickness+20;//=123;用于为Exi,Hyi分配数组  
const int		ii0=int(MX/2)-(Target_X+Interval_scatter_target); //total-scattered field boundary;ii0=10
const int		ii1=int(MX/2)+(Target_X+Interval_scatter_target);//ii1=56
const int		jj0=int(MY/2)-(Target_Y+Interval_scatter_target);//jj0=10　
const int		jj1=int(MY/2)+(Target_Y+Interval_scatter_target);//jj1=56　
	
//for PML in 2 dimension cartesian coordinates
const int		I1=PML_Thickness;
const int		I2=MX-PML_Thickness;　　
const int		J1=PML_Thickness;
const int		J2=MY-PML_Thickness;　　
const double	dte=dt/Epsilon;			//air  
const double	dtm=dt/Mu;                      
const double	me=Mu/Epsilon;                 

const double	R0=0.001;	//the reflect coefficent on the ObjectSpace-PML interface吸收边界处的反射系数
const int	index=3;	//variant index for Berenger's PML
const double	sigma_max=0.0-log(R0)*(index+1)/2.0*Epsilon*LightSpeed/PML_Thickness/__min(dx,dy);　　//?????????????????

double inline conduct(double L){
	if (L<=0) return 0.0;
	else return sigma_max*pow(L/PML_Thickness,index);   //返回sigma_e,PML电导率几何递增分布
}
double inline mconduct(double  L){
	return me*conduct(L+0.5);   //返回sigma_n,与sigma_e相差me=Mu/Epsilon,0.5，PML磁导率几何递增分布
}


//	the part for exciting source

double scp(int t){　　　　　　　　　　　　　　//高斯脉冲_源
	if(t<=0)  return 0.0;
	else
		return exp(0.0-(t-5.0*tao)*(t-5.0*tao)/(tao*tao/2.0));//tao=1e-10/dt
}

double scp0(int t){　　　　　　　　　　　　　//微分高斯脉冲_源
	if(0<=t&&t<=10*tao) 
		return -1000*(t-5.0*tao)*exp(0.0-(t-5.0*tao)*(t-5.0*tao)/(tao*tao/2.0));
	else return(0.0);
}

double sch(int t, const double aa=2.0){　　　//升余弦脉冲_源
	double enevop=1.0,w=2.0*pii*Frequency;
	int    T_max=(int)(0.5+aa/Frequency/dt);
	if (t<0) enevop=0.0;
	else if ((0<=t)&&(t<=T_max)) enevop=(0.5*(1.0-cos(w*t*dt*0.5/aa)));
		else enevop=1.0;		
	return(enevop*sin(w*t*dt));
}

double sins(int t){
	if(t<=0)  return 0.0;
	else      return(sin(2.0*pii*Frequency*t*dt));
}

//	the part for exciting source end



int inline _cdecl round(double x){

	int xx=(int)x;

	if(x>0.0){
		if(x-xx>=0.5)
			return (xx+1);
		else
			return xx;
	}
	else{
		if(fabs(x-xx)>=0.5)
			return(xx-1);
		else
			return xx;
	}
}