fdtd_3d_pml_hx.cpp

采用UMPL边界的FDTD编程

		for (j = nPML_y_1; j < ny_MIN_1_MIN_nPML_y_2; j++)
		{	
			ka = 0;
			for (k = nPML_z_1; k < nz_MIN_1_MIN_nPML_z_2; k++)
			{
				Bx_r = Bx_15[ia][ja][ka];

				Bx_15[ia][ja][ka] += ( (Ey[i][j][k+1] - Ey[i][j][k])*dtDIVdz -
							           (Ez[i][j+1][k] - Ez[i][j][k])*dtDIVdy );

				Hx[i][j][k] += inv_TwoOREp0*( K_Hx_c_2[ii]*Bx_15[ia][ja][ka] + 
							                  K_Hx_d_2[ii]*Bx_r )/Mat[Ind[i][j][k]][1];
				ka++;
			}
			ja++;
		}
		ia++;
		ii--;
	}
}

void CFDTD_3D_PML::calc_Hx_16()
{
	long  i, j, k, jj, ja, ka;
	double Bx_r;

	//16
	for (i = 0; i < nPML_x_1; i++)
	{	
		ja = 0;
		jj = nPML_y_2_MIN_1;
		for (j = ny_MIN_1_MIN_nPML_y_2; j < ny_MIN_1; j++)
		{	
			ka = 0;
			for (k = nPML_z_1; k < nz_MIN_1_MIN_nPML_z_2; k++)
			{
				Bx_r = Bx_16[i][ja][ka];

				Bx_16[i][ja][ka] = K_Bx_a_2[jj]*Bx_16[i][ja][ka] + 
					               K_Bx_b_2[jj]*( (Ey[i][j][k+1] - Ey[i][j][k])*inv_dz -
							                      (Ez[i][j+1][k] - Ez[i][j][k])*inv_dy);

				Hx[i][j][k] += inv_TwoOREp0*( K_Hx_c_1[i]*Bx_16[i][ja][ka] + 
							                  K_Hx_d_1[i]*Bx_r )/Mat[Ind[i][j][k]][1];
				ka++;
			}
			ja++;
			jj--;
		}
	}
}

void CFDTD_3D_PML::calc_Hx_17()
{
	long  i, j, k, jj;
	
	//17
	for (i = nPML_x_1; i < nx_MIN_nPML_x_2; i++)
	{	
		jj = nPML_y_2_MIN_1;
		for (j = ny_MIN_1_MIN_nPML_y_2; j < ny_MIN_1; j++)
		{	
			for (k = nPML_z_1; k < nz_MIN_1_MIN_nPML_z_2; k++)
			{
			   Hx[i][j][k] = K_Bx_a_2[jj]*Hx[i][j][k] + K_Bx_b_2[jj]/Mat[Ind[i][j][k]][1]*
					                 ( (Ey[i][j][k+1] - Ey[i][j][k])*inv_dzORmu_0 -
  						               (Ez[i][j+1][k] - Ez[i][j][k])*inv_dyORmu_0);
			}
			jj--;
		}
	}
}

void CFDTD_3D_PML::calc_Hx_18()
{
	long  i, j, k, ii, jj, ia, ja, ka;
	double Bx_r;

	//18
	ia = 0;
	ii = nPML_x_2_MIN_1;
	for (i = nx_MIN_nPML_x_2; i < nx; i++)
	{	
		ja = 0;
		jj = nPML_y_2_MIN_1;
		for (j = ny_MIN_1_MIN_nPML_y_2; j < ny_MIN_1; j++)
		{	
			ka = 0;
			for (k = nPML_z_1; k < nz_MIN_1_MIN_nPML_z_2; k++)
			{
				Bx_r = Bx_18[ia][ja][ka];

				Bx_18[ia][ja][ka] = K_Bx_a_2[jj]*Bx_18[ia][ja][ka] + 
					                K_Bx_b_2[jj]*( (Ey[i][j][k+1] - Ey[i][j][k])*inv_dz -
							                       (Ez[i][j+1][k] - Ez[i][j][k])*inv_dy);

				Hx[i][j][k] += inv_TwoOREp0*( K_Hx_c_2[ii]*Bx_18[ia][ja][ka] + 
							                  K_Hx_d_2[ii]*Bx_r )/Mat[Ind[i][j][k]][1];
				ka++;
			}
			ja++;
			jj--;
		}
		ia++;
		ii--;
	}
}

void CFDTD_3D_PML::calc_Hx_19()
{
	long  i, j, k, kk, ka;
	double Bx_r;

	//19
	for (i = 0; i < nPML_x_1; i++)
	{	
		for (j = 0; j < nPML_y_1; j++)
		{	
			ka = 0;
			kk = nPML_z_2_MIN_1;
			for (k = nz_MIN_1_MIN_nPML_z_2; k < nz_MIN_1; k++)
			{
				Bx_r = Bx_19[i][j][ka];

				Bx_19[i][j][ka] = K_Bx_a_1[j]*Bx_19[i][j][ka] + 
					              K_Bx_b_1[j]*( (Ey[i][j][k+1] - Ey[i][j][k])*inv_dz -
							                    (Ez[i][j+1][k] - Ez[i][j][k])*inv_dy);

				Hx[i][j][k] = K_Hx_a_2[kk]*Hx[i][j][k] + 
					          K_Hx_b_2[kk]*( K_Hx_c_1[i]*Bx_19[i][j][ka] + 
							                 K_Hx_d_1[i]*Bx_r )/Mat[Ind[i][j][k]][1];
				ka++;
				kk--;
			}
		}
	}
}

void CFDTD_3D_PML::calc_Hx_20()
{
	long  i, j, k, kk, ia, ka;
	double Bx_r;

	//20
	ia = 0;
	for (i = nPML_x_1; i < nx_MIN_nPML_x_2; i++)
	{	
		for (j = 0; j < nPML_y_1; j++)
		{	
			ka = 0;
			kk = nPML_z_2_MIN_1;
			for (k = nz_MIN_1_MIN_nPML_z_2; k < nz_MIN_1; k++)
			{
				Bx_r = Bx_20[ia][j][ka];

				Bx_20[ia][j][ka] = K_Bx_a_1[j]*Bx_20[ia][j][ka] + 
					               K_Bx_b_1[j]*( (Ey[i][j][k+1] - Ey[i][j][k])*inv_dz -
							                     (Ez[i][j+1][k] - Ez[i][j][k])*inv_dy);

				Hx[i][j][k] = K_Hx_a_2[kk]*Hx[i][j][k] + K_Hx_b_2[kk]*TwoOREp0DIVMu0*
					          ( Bx_20[ia][j][ka] - Bx_r )/Mat[Ind[i][j][k]][1];
				ka++;
				kk--;
			}
		}
		ia++;
	}
}

void CFDTD_3D_PML::calc_Hx_21()
{
	long  i, j, k, ii, kk, ia, ka;
	double Bx_r;

	//21
	ia = 0;
	ii = nPML_x_2_MIN_1;
	for (i = nx_MIN_nPML_x_2; i < nx; i++)
	{	
		for (j = 0; j < nPML_y_1; j++)
		{	
			ka = 0;
			kk = nPML_z_2_MIN_1;
			for (k = nz_MIN_1_MIN_nPML_z_2; k < nz_MIN_1; k++)
			{
				Bx_r = Bx_21[ia][j][ka];

				Bx_21[ia][j][ka] = K_Bx_a_1[j]*Bx_21[ia][j][ka] + 
					               K_Bx_b_1[j]*( (Ey[i][j][k+1] - Ey[i][j][k])*inv_dz -
							                     (Ez[i][j+1][k] - Ez[i][j][k])*inv_dy);

				Hx[i][j][k] = K_Hx_a_2[kk]*Hx[i][j][k] + 
					          K_Hx_b_2[kk]*( K_Hx_c_2[ii]*Bx_21[ia][j][ka] + 
							                 K_Hx_d_2[ii]*Bx_r )/Mat[Ind[i][j][k]][1];
				ka++;
				kk--;
			}
		}
		ia++;
		ii--;
	}
}

void CFDTD_3D_PML::calc_Hx_22()
{
	long  i, j, k, kk, ja, ka;
	double Bx_r;

	//22
	for (i = 0; i < nPML_x_1; i++)
	{	
		ja = 0;
		for (j = nPML_y_1; j < ny_MIN_1_MIN_nPML_y_2; j++)
		{	
			ka = 0;
			kk = nPML_z_2_MIN_1;
			for (k = nz_MIN_1_MIN_nPML_z_2; k < nz_MIN_1; k++)
			{
				Bx_r = Bx_22[i][ja][ka];

				Bx_22[i][ja][ka] += ( (Ey[i][j][k+1] - Ey[i][j][k])*dtDIVdz -
							          (Ez[i][j+1][k] - Ez[i][j][k])*dtDIVdy );

				Hx[i][j][k] = K_Hx_a_2[kk]*Hx[i][j][k] + 
					          K_Hx_b_2[kk]*( K_Hx_c_1[i]*Bx_22[i][ja][ka] + 
							                 K_Hx_d_1[i]*Bx_r )/Mat[Ind[i][j][k]][1];
				ka++;
				kk--;
			}
			ja++;
		}
	}
}

void CFDTD_3D_PML::calc_Hx_23()
{
	long  i, j, k, kk;
	
	//23
	for (i = nPML_x_1; i < nx_MIN_nPML_x_2; i++)
	{	
		for (j = nPML_y_1; j < ny_MIN_1_MIN_nPML_y_2; j++)
		{	
			kk = nPML_z_2_MIN_1;
			for (k = nz_MIN_1_MIN_nPML_z_2; k < nz_MIN_1; k++)
			{
			  Hx[i][j][k] = K_Hx_a_2[kk]*Hx[i][j][k] + K_Hx_b_2[kk]/Mat[Ind[i][j][k]][1]*
					       ( (Ey[i][j][k+1] - Ey[i][j][k])*TwoOREp0DIVMu0ORdtDIVdz -
					         (Ez[i][j+1][k] - Ez[i][j][k])*TwoOREp0DIVMu0ORdtDIVdy );
			  kk--;
			}
		}
	}
}

void CFDTD_3D_PML::calc_Hx_24()
{
	long  i, j, k, ii, kk, ia, ja, ka;
	double Bx_r;

	//24
	ia = 0;
	ii = nPML_x_2_MIN_1;
	for (i = nx_MIN_nPML_x_2; i < nx; i++)
	{	
		ja = 0;
		for (j = nPML_y_1; j < ny_MIN_1_MIN_nPML_y_2; j++)
		{	
			ka = 0;
			kk = nPML_z_2_MIN_1;
			for (k = nz_MIN_1_MIN_nPML_z_2; k < nz_MIN_1; k++)
			{
				Bx_r = Bx_24[ia][ja][ka];

				Bx_24[ia][ja][ka] += ( (Ey[i][j][k+1] - Ey[i][j][k])*dtDIVdz -
							           (Ez[i][j+1][k] - Ez[i][j][k])*dtDIVdy );

				Hx[i][j][k] = K_Hx_a_2[kk]*Hx[i][j][k] + 
					          K_Hx_b_2[kk]*( K_Hx_c_2[ii]*Bx_24[ia][ja][ka] + 
							                K_Hx_d_2[ii]*Bx_r )/Mat[Ind[i][j][k]][1];
				ka++;
				kk--;
			}
			ja++;
		}
		ia++;
		ii--;
	}
}

void CFDTD_3D_PML::calc_Hx_25()
{
	long  i, j, k, jj, kk, ja, ka;
	double Bx_r;

	//25
	for (i = 0; i < nPML_x_1; i++)
	{	
		ja = 0;
		jj = nPML_y_2_MIN_1;
		for (j = ny_MIN_1_MIN_nPML_y_2; j < ny_MIN_1; j++)
		{	
			ka = 0;
			kk = nPML_z_2_MIN_1;
			for (k = nz_MIN_1_MIN_nPML_z_2; k < nz_MIN_1; k++)
			{
				Bx_r = Bx_25[i][ja][ka];

				Bx_25[i][ja][ka] = K_Bx_a_2[jj]*Bx_25[i][ja][ka] + 
					               K_Bx_b_2[jj]*( (Ey[i][j][k+1] - Ey[i][j][k])*inv_dz -
							                      (Ez[i][j+1][k] - Ez[i][j][k])*inv_dy);

				Hx[i][j][k] = K_Hx_a_2[kk]*Hx[i][j][k] + 
					          K_Hx_b_2[kk]*( K_Hx_c_1[i]*Bx_25[i][ja][ka] + 
							                 K_Hx_d_1[i]*Bx_r )/Mat[Ind[i][j][k]][1];
				ka++;
				kk--;
			}
			ja++;
			jj--;
		}
	}
}

void CFDTD_3D_PML::calc_Hx_26()
{
	long  i, j, k, jj, kk, ia, ja, ka;
	double Bx_r;

	//26
	ia = 0;
	for (i = nPML_x_1; i < nx_MIN_nPML_x_2; i++)
	{	
		ja = 0;
		jj = nPML_y_2_MIN_1;
		for (j = ny_MIN_1_MIN_nPML_y_2; j < ny_MIN_1; j++)
		{	
			ka = 0;
			kk = nPML_z_2_MIN_1;
			for (k = nz_MIN_1_MIN_nPML_z_2; k < nz_MIN_1; k++)
			{
				Bx_r = Bx_26[ia][ja][ka];

				Bx_26[ia][ja][ka] = K_Bx_a_2[jj]*Bx_26[ia][ja][ka] + 
					                K_Bx_b_2[jj]*( (Ey[i][j][k+1] - Ey[i][j][k])*inv_dz -
							                       (Ez[i][j+1][k] - Ez[i][j][k])*inv_dy);

				Hx[i][j][k] = K_Hx_a_2[kk]*Hx[i][j][k] + K_Hx_b_2[kk]*TwoOREp0DIVMu0*
					           ( Bx_26[ia][ja][ka] - Bx_r )/Mat[Ind[i][j][k]][1];
				ka++;
				kk--;
			}
			ja++;
			jj--;
		}
		ia++;
	}
}

void CFDTD_3D_PML::calc_Hx_27()
{
	long  i, j, k, ii, jj, kk, ia, ja, ka;
	double Bx_r;

	//27
	ia = 0;
	ii = nPML_x_2_MIN_1;
	for (i = nx_MIN_nPML_x_2; i < nx; i++)
	{	
		ja = 0;
		jj = nPML_y_2_MIN_1;
		for (j = ny_MIN_1_MIN_nPML_y_2; j < ny_MIN_1; j++)
		{	
			ka = 0;
			kk = nPML_z_2_MIN_1;
			for (k = nz_MIN_1_MIN_nPML_z_2; k < nz_MIN_1; k++)
			{
				

				Bx_r = Bx_27[ia][ja][ka];

				Bx_27[ia][ja][ka] = K_Bx_a_2[jj]*Bx_27[ia][ja][ka] + 
					                K_Bx_b_2[jj]*( (Ey[i][j][k+1] - Ey[i][j][k])*inv_dz -
							                       (Ez[i][j+1][k] - Ez[i][j][k])*inv_dy);

				Hx[i][j][k] = K_Hx_a_2[kk]*Hx[i][j][k] + 
					          K_Hx_b_2[kk]*( K_Hx_c_2[ii]*Bx_27[ia][ja][ka] + 
							                 K_Hx_d_2[ii]*Bx_r )/Mat[Ind[i][j][k]][1];
				ka++;
				kk--;
			}
			ja++;
			jj--;
		}
		ia++;
		ii--;
	}
}