fdtd_xypbc_zpml.h

利用c++语言写的三维FDTD

#pragma once

#include "stdafx.h"
#include "FDTD_1D_ExHy.h"

class CFDTD_xyPBC_zPML
{
private:
	double pi, TwoORpi, eps_0, mu_0;

	int nr_threads;

	long *** Ind; //indicates the material type in each FDTD cell
	double ** Mat; //contains the material properties  eps_r sigma

	long **Coord_ptSource; //the location of the point sources
	double **Par_ptSource; //the parameters of the point sources

	//The field components
	double ***Ex, ***Ey, ***Fz, ***Ez, ***Hx, ***Hy, ***Gz, ***Hz;

	//Coefficients containing the PML boundary parameters
	double *K_E1_a, *K_E1_b, *K_E3_a, *K_E3_b, *K_E4_a, *K_E4_b, *K_E6_a, *K_E6_b;
    	
	//dimensions of the computational space
    long nx, ny, nz;
	//dimensions of the PML region
	long n_z_l, n_z_h;

	//number of materials present in the investigated geometry
	long n_Mat;
	
	//number of point sources present in the computational domain
	long n_ptSource;
	long source_type;
	long jel_plane_wave;

	long n_TS_z, n_1D; //the z position of the total field scatterd field interface
	
	CFDTD_1D_ExHy fdtd_1D;
	
	double *Ex_1D, *Hy_1D;

	//time increment
	double dt;
	//elementary cell sizes
	double dx, dy, dz; 
	
	double inv_dx, inv_dy, inv_dz;

	//the coordinates of the followed field components
	long **Ind_Foll; 
	 //number of the followed points and the number of iterations
	long length_Ind_Foll, n_tot;
	//stores the followed field components
	double **Hx_Foll, **Hy_Foll, **Hz_Foll, **Ex_Foll, **Ey_Foll, **Ez_Foll;
	
	long ***Mask_Wigner;

	//stores the averaged field components
	int aver_field_volume;
	double *Ex_Foll_av, *Ey_Foll_av, *Ez_Foll_av;
	double *Hx_Foll_av, *Hy_Foll_av, *Hz_Foll_av;
	long nx_a_av, nx_b_av, ny_a_av, ny_b_av, nz_a_av, nz_b_av;
	long Vol_av;
	
	/////////////////////////
	//Fourier transform in a subvolume
	/////////////////////////
	int jel_fourier_transf_vol;
	long nx_a_f, nx_b_f, ny_a_f, ny_b_f, nz_a_f, nz_b_f;
	double frec_1, frec_2;
	long n_frec, n_f_Points, start_fourier;
	double **Ex_fourier_TR_re, **Ex_fourier_TR_im; 
	double **Ey_fourier_TR_re, **Ey_fourier_TR_im; 
	double **Ez_fourier_TR_re, **Ez_fourier_TR_im; 
	double **Hx_fourier_TR_re, **Hx_fourier_TR_im; 
	double **Hy_fourier_TR_re, **Hy_fourier_TR_im; 
	double **Hz_fourier_TR_re, **Hz_fourier_TR_im; 
	double *fourier_omega;

	//the Bragg diffraction
	int jel_aver_Bragg_WignerUC, nz_Bragg;
	double kx, ky;
	double **Gxy;
	long n_Gxy;
	double **Ex_Br_av_re, **Ey_Br_av_re, **Ez_Br_av_re;
	double **Hx_Br_av_re, **Hy_Br_av_re, **Hz_Br_av_re;
	double **Ex_Br_av_im, **Ey_Br_av_im, **Ez_Br_av_im;
	double **Hx_Br_av_im, **Hy_Br_av_im, **Hz_Br_av_im;
	long Surf_Bragg;
	
	//to save the field components
	long n_x_a, n_x_b, n_y_a, n_y_b, n_z_a, n_z_b; //identify the saved volume
	long nx_yz, ny_xz, nz_xy;
	char *path_name_Ex;
	char *path_name_Ey;
	char *path_name_Ez;
	char *path_name_Hx;
	char *path_name_Hy;
	char *path_name_Hz;
	char *path_name_Ex_1D;
	int jel_Save_Slice;
	long saveFROMinst, saveTOinst;
	int nr_save;
	double NORMA;

	#ifndef run_enviroment //IBM AIX
		char *Path_Load_Workspace_Data, *Path_Save_Workspace_Data;

		//to measure the ellapsed time UNIX IBM_AIX
		double time_start, el_time, limit_time;
		struct timeval tv;
		struct timezone tz;
		int load_workspace_data;
	#endif

public:
	CFDTD_xyPBC_zPML(void);
	~CFDTD_xyPBC_zPML(void);
	
	void Init_nr_THR(int nr_Threads);

	bool Init(long ***&Index, long n_x, long n_y, long n_z, double **&Mater, long nMat, 
			  long n_PML, double d_t, double d_x, double d_y, double d_z);
	//Index - describe the geometry, each element indicating the material type in 
	//        the computational grid 
	//nx, ny, nz - the dimensions of the computational space
	//Mater - contains the permittivity and conductivity. Dimension [n_Mat][2] 
	//n_PML - indicate the dimension of the PML
	
	void Init_PML_Par(double eps_r_a, double mu_r_a,
		              double eps_r_b, double mu_r_b);
	void Init_plWave_Gauss(double E0, double t0, double tw, long n_z_TS);
	void Init_plWave_Sin(double E0, double omega, long n_z_TS);
	void Init_plWave_GaussSin(double E0, double omega, double t0, double tw, long n_z_TS);

	void Init_ptSource(long **Coord, long nCoord,  double **Par_Excit, long s_type);

	int fdtd_marching(long n_max_steep);
	
	//functions to save output data
	bool Init_Followed(long **Ind_Followed, long n, long n_t);
	int Init_Aver_xyz(long Nx_a_av, long Nx_b_av, long Ny_a_av, 
		              long Ny_b_av, long Nz_a_av, long Nz_b_av, 
					  long N_iter);

	bool Init_Save_FieldSlices(long nxa, long nxb, long nya, long nyb, 
						       long nza, long nzb, long nZ_xy, long nY_xz, 
							   long nX_yz, long t1, long t2, long nr_Save, 
							   char *path);

	void Get_Data_Followed(double **&h_x, double **&h_y, double **&h_z, 
						   double **&e_x, double **&e_y, double **&e_z);

	void Get_Data_Aver(double *&hx_av, double *&hy_av, double *&hz_av, 
					   double *&ex_av, double *&ey_av, double *&ez_av);


	int aver_xyz_WignerC(double ***X, long nx_a_av, long nx_b_av, long ny_a_av, long ny_b_av,
			             long nz_a_av, long nz_b_av, double *X_av, double Vol);

	int Init_aver_Bragg_WignerC(long nZ_av, double K_x, double K_y, 
		                        double **&Gxy, long n_Gxy, long N_iter);

	int aver_Bragg_WignerC(double ***X, double kx, double ky, 
						   double Gx, double Gy, double *X_avBr_re, 
						   double *X_avBr_im, double Surf);
	
	void Get_Data_Bragg_Aver(double **&ex_avBr_re, double **&ex_avBr_im, double **&ey_avBr_re, 
		                     double **&ey_avBr_im, double **&ez_avBr_re, double **&ez_avBr_im,
							 double **&hx_avBr_re, double **&hx_avBr_im, double **&hy_avBr_re, 
							 double **&hy_avBr_im, double **&hz_avBr_re, double **&hz_avBr_im);

	void Init_Mask_Wigner(long ***&mask, long nx_W, long ny_W, long nz_W);

	int Init_Fourier_Transf(long Nx_a_f, long Nx_b_f, long Ny_a_f, long Ny_b_f,
							long Nz_a_f, long Nz_b_f, double Frec_1, double Frec_2,
	                        long N_frec, long Start_fourier);

	int fourier_transf_vol(double ***U, long nx_a_f, long nx_b_f, long ny_a_f, long ny_b_f, 
				           long nz_a_f, long nz_b_f, double **FrT_re, double **FrT_im, 
				           double *FrT_om, long n_FrT_om, double time);

	void Get_Data_Fourier(double **&ex_fourier_TR_re, double **&ex_fourier_TR_im,
						  double **&ey_fourier_TR_re, double **&ey_fourier_TR_im,
						  double **&ez_fourier_TR_re, double **&ez_fourier_TR_im,
						  double **&hx_fourier_TR_re, double **&hx_fourier_TR_im,
						  double **&hy_fourier_TR_re, double **&hy_fourier_TR_im,
						  double **&hz_fourier_TR_re, double **&hz_fourier_TR_im,
						  double *&omega);

	#ifndef run_enviroment //IBM AIX
		int Init_Mesure_Ellapsed_Time(double TimeStart, double LimitTime, int LoadWorkspaceData, 
									  char *path_load, char *path_save);
		int load_FE_BH(char *Path);
		int save_FE_BH(long iter, double time, char *Path);
	#endif

	//free up the allocated memory	
	void Free_Mem();
	
private:
	//the field components from Maxwell equations
	void calc_Ex();
	void calc_Ey();
	void calc_Ez();
	void calc_Hx();
	void calc_Hy();
	void calc_Hz();

	//Periodic BC
	void PBC_Ex();
	void PBC_Ey();
	void PBC_Ez();

	void PBC_Hx();
	void PBC_Hy();
	void PBC_Hz();

	//Functions referring to the excitation
	//total field scattered field interface in z direction
	void TS_Ex();
	void TS_Hy();
	//Point Source
	void ptSource_J(double ***&xx, double **&Par_ptS, double time);

	//the coordinates of the followed field components
	void Set_Data_Followed(long n_t);
	int Save_FieldSlices(long iter);
	int Save_Field(long iter);
};