⭐ 欢迎来到虫虫下载站！ | 📦 资源下载 📁 资源专辑 ℹ️ 关于我们
⭐ 虫虫下载站
📤 上传资源

📄 xiangyang27.cpp

📁 3D-FDTD垂直双阵子程序
💻 CPP
📖 第 1 页 / 共 4 页
字号:
🔍
1234 
	 {  for(j=0;j<=180;j++)
	      fprintf(out,"%f\n",sqrt(fxhs[j].shi*fxhs[j].shi+fxhs[j].xu*fxhs[j].xu));
     }

end:  return(0);
 }







void fuzhi()
{int loop1,loop2,loop3;
    for(loop1=0;loop1<Nx;loop1++)
	{ for(loop2=0;loop2<Ny+1;loop2++)
	   {for(loop3=0;loop3<Nz+1;loop3++)
	{      Exy[loop1][loop2][loop3]=0;
	       Exz[loop1][loop2][loop3]=0;
	}
	}
    }   /*给x向电场赋初值*/

    for(loop1=0;loop1<Nx+1;loop1++)
	{ for(loop2=0;loop2<Ny;loop2++)
	   {for(loop3=0;loop3<Nz+1;loop3++)
	{     Eyx[loop1][loop2][loop3]=0;
	      Eyz[loop1][loop2][loop3]=0;
	}
	}
    }   /*给y向电场赋初值*/

   for(loop1=0;loop1<Nx+1;loop1++)
	{ for(loop2=0;loop2<Ny+1;loop2++)
	   {for(loop3=0;loop3<Nz;loop3++)
     
   {       Ezx[loop1][loop2][loop3]=0;
	       Ezy[loop1][loop2][loop3]=0;
   }
       
	}
    }   /*给z向电场赋初值,并加激励源*/

      for(loop1=0;loop1<Nx+1;loop1++)
	{ for(loop2=0;loop2<Ny;loop2++)
	   {for(loop3=0;loop3<Nz;loop3++)
	{Hxy[loop1][loop2][loop3]=0;
     Hxz[loop1][loop2][loop3]=0;

	}
	}
    }   /*给x向磁场赋初值*/

     for(loop1=0;loop1<Nx;loop1++)
   { for(loop2=0;loop2<Ny+1;loop2++)
	   {for(loop3=0;loop3<Nz;loop3++)
	{Hyx[loop1][loop2][loop3]=0;
     Hyz[loop1][loop2][loop3]=0;
	}
	}
    }   /*给y向磁场赋初值*/

      for(loop1=0;loop1<Nx;loop1++)
	{ for(loop2=0;loop2<Ny;loop2++)
	   {for(loop3=0;loop3<Nz+1;loop3++)
	{  Hzx[loop1][loop2][loop3]=0;
       Hzy[loop1][loop2][loop3]=0;
	}
	}
    }   /*给z向磁场赋初值*/
    

	  for(loop1=0;loop1<NN/2;loop1++)
	  {
	    pp1[loop1].shi=0;pp1[loop1].xu=0;
		 pp2[loop1].shi=0;pp2[loop1].xu=0;
		  zk[loop1].shi=0;zk[loop1].xu=0;
	  }

	   for(loop1=0;loop1<=180;loop1++)
	   {fxhs[loop1].shi=0,fxhs[loop1].xu=0;}

		for(loop1=0;loop1<=180;loop1++)
		{si[loop1]=sin(loop1/180.0*Pi);co[loop1]=cos(loop1/180.0*Pi);}

		 for(loop1=0;loop1<Nx-2*NNN-2*DDD;loop1++)
			for(loop2=0;loop2<Nx-2*NNN-2*DDD+1;loop2++)
			{
			   Jxu[loop1][loop2].zhengfu=0; Jxu[loop1][loop2].xiangwei=0;
			   Jxd[loop1][loop2].zhengfu=0; Jxd[loop1][loop2].xiangwei=0;
			   Jxl[loop1][loop2].zhengfu=0; Jxl[loop1][loop2].xiangwei=0;
			   Jxr[loop1][loop2].zhengfu=0; Jxr[loop1][loop2].xiangwei=0;
			   Jzl[loop1][loop2].zhengfu=0; Jzl[loop1][loop2].xiangwei=0;
			   Jzr[loop1][loop2].zhengfu=0; Jzr[loop1][loop2].xiangwei=0;
			   Jmxu[loop1][loop2].zhengfu=0; Jmxu[loop1][loop2].xiangwei=0;
			   Jmxd[loop1][loop2].zhengfu=0; Jmxd[loop1][loop2].xiangwei=0;
			   Jmxl[loop1][loop2].zhengfu=0; Jmxl[loop1][loop2].xiangwei=0;
			   Jmxr[loop1][loop2].zhengfu=0; Jmxr[loop1][loop2].xiangwei=0;
			   Jmzl[loop1][loop2].zhengfu=0; Jmzl[loop1][loop2].xiangwei=0;
			   Jmzr[loop1][loop2].zhengfu=0; Jmzr[loop1][loop2].xiangwei=0;
			}

	    for(loop1=0;loop1<Nx-2*NNN-2*DDD+1;loop1++)
			for(loop2=0;loop2<Nx-2*NNN-2*DDD;loop2++)
			{
			   Jyu[loop1][loop2].zhengfu=0; Jyu[loop1][loop2].xiangwei=0;
			   Jyd[loop1][loop2].zhengfu=0; Jyd[loop1][loop2].xiangwei=0;
			   Jyf[loop1][loop2].zhengfu=0; Jyf[loop1][loop2].xiangwei=0;
			   Jyb[loop1][loop2].zhengfu=0; Jyb[loop1][loop2].xiangwei=0;
			   Jzf[loop1][loop2].zhengfu=0; Jzf[loop1][loop2].xiangwei=0;
			   Jzb[loop1][loop2].zhengfu=0; Jzb[loop1][loop2].xiangwei=0;
			   Jmyu[loop1][loop2].zhengfu=0; Jmyu[loop1][loop2].xiangwei=0;
			   Jmyd[loop1][loop2].zhengfu=0; Jmyd[loop1][loop2].xiangwei=0;
			   Jmyf[loop1][loop2].zhengfu=0; Jmyf[loop1][loop2].xiangwei=0;
			   Jmyb[loop1][loop2].zhengfu=0; Jmyb[loop1][loop2].xiangwei=0;
			   Jmzf[loop1][loop2].zhengfu=0; Jmzf[loop1][loop2].xiangwei=0;
			   Jmzb[loop1][loop2].zhengfu=0; Jmzb[loop1][loop2].xiangwei=0;
			}
}



/***************************************************************************************************************************************/
 /*x向电场*/
   int chuanbo1(int M)  /*设置波动源，向前波动1个周期，并实现前时刻和后一时刻值的传递。M为时间点，需要返回的参数以便进行调试*/
   {int loop1,loop2,loop3;
      for(loop1=0;loop1<Nx;loop1++)
	  {for(loop2=1;loop2<Ny;loop2++)
	  {for(loop3=1;loop3<Nz;loop3++)
      //if(loop1==N-NNN-5&&loop2==N-NNN-5&&loop3==N-NNN-5)
/*		if(loop1==Nx/2-DD/2+kx&&loop2==Ny/2-DD/2+ky&&loop3>=NNN+5&&loop3<NNN+5+5)   //前一个5是预留网格,后一个5是介质厚度
	   {
		   Exy[loop1][loop2][loop3]=((2*e*er-sig*dt)/(2*e*er+sig*dt))*Exy[loop1][loop2][loop3]+(2*dt)/(2*e*er+sig*dt)*((Hzx[loop1][loop2][loop3]+Hzy[loop1][loop2][loop3]-Hzx[loop1][loop2-1][loop3]-Hzy[loop1][loop2-1][loop3])/dy)-sin(2*Pi*M*dt*f);
	       Exz[loop1][loop2][loop3]=((2*e*er-sig*dt)/(2*e*er+sig*dt))*Exz[loop1][loop2][loop3]-(2*dt)/(2*e*er+sig*dt)*((Hyx[loop1][loop2][loop3]+Hyz[loop1][loop2][loop3]-Hyz[loop1][loop2][loop3-1]-Hyx[loop1][loop2][loop3-1])/dz);
	   }*/
	  {   if(loop1>NNN-1&&loop1<Nx-NNN&&loop2>NNN-1&&loop2<Ny-NNN&&loop3>NNN-1&&loop3<Nz-NNN)
		
		  if(loop1>=0&&loop1<=NNN-1)
		 {  sgmx=sgmmax*(NNN-loop1-1)/NNN*(NNN-loop1-1)/NNN;
		    sgmmx=sgmx*142122.3034;
		 }
		  else if(loop1>=Nx-NNN&&loop1<Nx)
		 {  sgmx=sgmmax*(NNN+loop1-Nx)/NNN*(NNN+loop1-Nx)/NNN;
		    sgmmx=sgmx*142122.3034;
		 }
		  else
		  {sgmx=0;sgmmx=0;}


           if(loop2>=0&&loop2<=NNN-1)
		 {  sgmy=sgmmax*(NNN-loop2-1)/NNN*(NNN-loop2-1)/NNN;
		    sgmmy=sgmy*142122.3034;
		 }
		  else if(loop2>=Ny-NNN&&loop2<Ny)
		 {  sgmy=sgmmax*(NNN+loop2-Ny)/NNN*(NNN+loop2-Ny)/NNN;
		    sgmmy=sgmy*142122.3034;
		 }
		  else
		  {sgmy=0;sgmmy=0;}


		   if(loop3>=0&&loop3<=NNN-1)
		 {  sgmz=sgmmax*(NNN-loop3-1)/NNN*(NNN-loop3-1)/NNN;
		    sgmmz=sgmz*142122.3034;
		 }
		  else if(loop3>=Nz-NNN&&loop3<Nz)
		 {  sgmz=sgmmax*(NNN+loop3-Nz)/NNN*(NNN+loop3-Nz)/NNN;
		    sgmmz=sgmz*142122.3034;
		 }
		  else
		  {sgmz=0;sgmmz=0;}


   
      Exy[loop1][loop2][loop3]=((2*e-sgmy*dt)/(2*e+sgmy*dt))*Exy[loop1][loop2][loop3]+(2*dt)/(2*e+sgmy*dt)*((Hzx[loop1][loop2][loop3]+Hzy[loop1][loop2][loop3]-Hzx[loop1][loop2-1][loop3]-Hzy[loop1][loop2-1][loop3])/dy);
	  Exz[loop1][loop2][loop3]=((2*e-sgmz*dt)/(2*e+sgmz*dt))*Exz[loop1][loop2][loop3]-(2*dt)/(2*e+sgmz*dt)*((Hyx[loop1][loop2][loop3]+Hyz[loop1][loop2][loop3]-Hyz[loop1][loop2][loop3-1]-Hyx[loop1][loop2][loop3-1])/dz);
	  if(loop1>NNN-1&&loop1<Nx-NNN&&loop2>NNN-1&&loop2<Ny-NNN&&loop3>NNN-1&&loop3<Nz-NNN)
	  { Exy[loop1][loop2][loop3]+= Exz[loop1][loop2][loop3];
	    Exz[loop1][loop2][loop3]=0;
	  }
	  e=ee;
	   }
     }
	  
	          
	  }
	 
	 
/***********************************************************************************************************************************************************************************************************************************/
    /*y向电场*/
	  for(loop1=1;loop1<Nx;loop1++)
	  {for(loop2=0;loop2<Ny;loop2++)
	  {for(loop3=1;loop3<Nz;loop3++)
  /*    if(loop1==N/2&&loop2==N/2&&loop3==N/2)
	   {  Eyz[loop1][loop2][loop3]=((2*e-sgmz*dt)/(2*e+sgmz*dt))*Eyz[loop1][loop2][loop3]+(2*dt)/(2*e+sgmz*dt)*((Hxy[loop1][loop2][loop3]+Hxz[loop1][loop2][loop3]-Hxy[loop1][loop2][loop3-1]-Hxz[loop1][loop2][loop3-1])/dz)-sin(2*Pi*f*dt*(M-0.5));
		  Eyx[loop1][loop2][loop3]=((2*e-sgmx*dt)/(2*e+sgmx*dt))*Eyx[loop1][loop2][loop3]-(2*dt)/(2*e+sgmx*dt)*((Hzx[loop1][loop2][loop3]+Hzy[loop1][loop2][loop3]-Hzx[loop1-1][loop2][loop3]-Hzy[loop1-1][loop2][loop3])/dx);
	   }
 
	   else
	   */
	  {   	    
		  if(loop1>=0&&loop1<=NNN-1)
		 {  sgmx=sgmmax*(NNN-loop1-1)/NNN*(NNN-loop1-1)/NNN;
		    sgmmx=sgmx*142122.3034;
		 }
		  else if(loop1>=Nx-NNN&&loop1<Nx)
		 {  sgmx=sgmmax*(NNN+loop1-Nx)/NNN*(NNN+loop1-Nx)/NNN;
		    sgmmx=sgmx*142122.3034;
		 }
		  else
		  {sgmx=0;sgmmx=0;}


           if(loop2>=0&&loop2<=NNN-1)
		 {  sgmy=sgmmax*(NNN-loop2-1)/NNN*(NNN-loop2-1)/NNN;
		    sgmmy=sgmy*142122.3034;
		 }
		  else if(loop2>=Ny-NNN&&loop2<Ny)
		 {  sgmy=sgmmax*(NNN+loop2-Ny)/NNN*(NNN+loop2-Ny)/NNN;
		    sgmmy=sgmy*142122.3034;
		 }
		  else
		  {sgmy=0;sgmmy=0;}


		   if(loop3>=0&&loop3<=NNN-1)
		 {  sgmz=sgmmax*(NNN-loop3-1)/NNN*(NNN-loop3-1)/NNN;
		    sgmmz=sgmz*142122.3034;
		 }
		  else if(loop3>=Nz-NNN&&loop3<Nz)
		 {  sgmz=sgmmax*(NNN+loop3-Nz)/NNN*(NNN+loop3-Nz)/NNN;
		    sgmmz=sgmz*142122.3034;
		 }
		  else
		  {sgmz=0;sgmmz=0;}



        Eyz[loop1][loop2][loop3]=((2*e-sgmz*dt)/(2*e+sgmz*dt))*Eyz[loop1][loop2][loop3]+(2*dt)/(2*e+sgmz*dt)*((Hxy[loop1][loop2][loop3]+Hxz[loop1][loop2][loop3]-Hxy[loop1][loop2][loop3-1]-Hxz[loop1][loop2][loop3-1])/dz);
	    Eyx[loop1][loop2][loop3]=((2*e-sgmx*dt)/(2*e+sgmx*dt))*Eyx[loop1][loop2][loop3]-(2*dt)/(2*e+sgmx*dt)*((Hzx[loop1][loop2][loop3]+Hzy[loop1][loop2][loop3]-Hzx[loop1-1][loop2][loop3]-Hzy[loop1-1][loop2][loop3])/dx);
  
		   if(loop1>NNN-1&&loop1<Nx-NNN&&loop2>NNN-1&&loop2<Ny-NNN&&loop3>NNN-1&&loop3<Nz-NNN)
		   {Eyz[loop1][loop2][loop3]+= Eyx[loop1][loop2][loop3];
		    Eyx[loop1][loop2][loop3]=0;}
		   e=ee;
		
	  }
	  }
	  }
	 
/****************************************************************************************************************************************************************************************************************************************/      
 /*z向电场*/
	  for(loop1=0;loop1<Nx+1;loop1++)
	  {for(loop2=0;loop2<Ny+1;loop2++)
	  {for(loop3=0;loop3<Nz;loop3++)
		if(loop1==Nx/2&&loop2==Ny/2&&loop3==Nz/2+D/2+DD/4)   //前一个5是预留网格,后一个5是介质厚度
	   {
		  Ezx[loop1][loop2][loop3]=((2*e*er-sig*dt)/(2*e*er+sig*dt))*Ezx[loop1][loop2][loop3]+(2*dt)/(2*e*er+sig*dt)*((Hyz[loop1][loop2][loop3]+Hyx[loop1][loop2][loop3]-Hyz[loop1-1][loop2][loop3]-Hyx[loop1-1][loop2][loop3])/dx)+sin(2*Pi*f*(M*dt-0.5*dt)+1.5*Pi);
	//	  Ezx[loop1][loop2][loop3]=((2*e*er-sig*dt)/(2*e*er+sig*dt))*Ezx[loop1][loop2][loop3]+(2*dt)/(2*e*er+sig*dt)*((Hyz[loop1][loop2][loop3]+Hyx[loop1][loop2][loop3]-Hyz[loop1-1][loop2][loop3]-Hyx[loop1-1][loop2][loop3])/dx)-(sin(2*Pi*(M-0.5)*f*dt)+sin(2*Pi*(M-0.5)*(f+20000000)*dt)+sin(2*Pi*(M-0.5)*(f+40000000)*dt)+sin(2*Pi*(M-0.5)*(f+60000000)*dt)+sin(2*Pi*(M-0.5)*(f+80000000)*dt)+sin(2*Pi*(M-0.5)*(f+100000000)*dt)+sin(2*Pi*(M-0.5)*(f+120000000)*dt)+sin(2*Pi*(M-0.5)*(f+140000000)*dt)+sin(2*Pi*(M-0.5)*(f+160000000)*dt)+sin(2*Pi*(M-0.5)*(f+180000000)*dt)+sin(2*Pi*(M-0.5)*(f+200000000)*dt)+sin(2*Pi*(M-0.5)*(f+220000000)*dt)+sin(2*Pi*(M-0.5)*(f+240000000)*dt)+sin(2*Pi*(M-0.5)*(f+260000000)*dt)+sin(2*Pi*(M-0.5)*(f+280000000)*dt)+sin(2*Pi*(M-0.5)*(f+300000000)*dt));
		  Ezy[loop1][loop2][loop3]=((2*e*er-sig*dt)/(2*e*er+sig*dt))*Ezy[loop1][loop2][loop3]-(2*dt)/(2*e*er+sig*dt)*((Hxy[loop1][loop2][loop3]+Hxz[loop1][loop2][loop3]-Hxy[loop1][loop2-1][loop3]-Hxz[loop1][loop2-1][loop3])/dy);
		}
		else if(loop1==Nx/2&&loop2==Ny/2&&loop3==Nz/2-D/2-DD/4)   //前一个5是预留网格,后一个5是介质厚度
	   {
		  Ezx[loop1][loop2][loop3]=((2*e*er-sig*dt)/(2*e*er+sig*dt))*Ezx[loop1][loop2][loop3]+(2*dt)/(2*e*er+sig*dt)*((Hyz[loop1][loop2][loop3]+Hyx[loop1][loop2][loop3]-Hyz[loop1-1][loop2][loop3]-Hyx[loop1-1][loop2][loop3])/dx)+sin(2*Pi*f*(M*dt-0.5*dt));
	//	  Ezx[loop1][loop2][loop3]=((2*e*er-sig*dt)/(2*e*er+sig*dt))*Ezx[loop1][loop2][loop3]+(2*dt)/(2*e*er+sig*dt)*((Hyz[loop1][loop2][loop3]+Hyx[loop1][loop2][loop3]-Hyz[loop1-1][loop2][loop3]-Hyx[loop1-1][loop2][loop3])/dx)-(sin(2*Pi*(M-0.5)*f*dt)+sin(2*Pi*(M-0.5)*(f+20000000)*dt)+sin(2*Pi*(M-0.5)*(f+40000000)*dt)+sin(2*Pi*(M-0.5)*(f+60000000)*dt)+sin(2*Pi*(M-0.5)*(f+80000000)*dt)+sin(2*Pi*(M-0.5)*(f+100000000)*dt)+sin(2*Pi*(M-0.5)*(f+120000000)*dt)+sin(2*Pi*(M-0.5)*(f+140000000)*dt)+sin(2*Pi*(M-0.5)*(f+160000000)*dt)+sin(2*Pi*(M-0.5)*(f+180000000)*dt)+sin(2*Pi*(M-0.5)*(f+200000000)*dt)+sin(2*Pi*(M-0.5)*(f+220000000)*dt)+sin(2*Pi*(M-0.5)*(f+240000000)*dt)+sin(2*Pi*(M-0.5)*(f+260000000)*dt)+sin(2*Pi*(M-0.5)*(f+280000000)*dt)+sin(2*Pi*(M-0.5)*(f+300000000)*dt));
		  Ezy[loop1][loop2][loop3]=((2*e*er-sig*dt)/(2*e*er+sig*dt))*Ezy[loop1][loop2][loop3]-(2*dt)/(2*e*er+sig*dt)*((Hxy[loop1][loop2][loop3]+Hxz[loop1][loop2][loop3]-Hxy[loop1][loop2-1][loop3]-Hxz[loop1][loop2-1][loop3])/dy);
		}
		else
	  {   
		  if(loop1>=0&&loop1<=NNN-1)
		 {  sgmx=sgmmax*(NNN-loop1-1)/NNN*(NNN-loop1-1)/NNN;
		    sgmmx=sgmx*142122.3034;
		 }
		  else if(loop1>=Nx-NNN&&loop1<Nx)
		 {  sgmx=sgmmax*(NNN+loop1-Nx)/NNN*(NNN+loop1-Nx)/NNN;
		    sgmmx=sgmx*142122.3034;
		 }
		  else
		  {sgmx=0;sgmmx=0;}


           if(loop2>=0&&loop2<=NNN-1)
		 {  sgmy=sgmmax*(NNN-loop2-1)/NNN*(NNN-loop2-1)/NNN;
		    sgmmy=sgmy*142122.3034;
		 }
		  else if(loop2>=Ny-NNN&&loop2<Ny)
		 {  sgmy=sgmmax*(NNN+loop2-Ny)/NNN*(NNN+loop2-Ny)/NNN;
		    sgmmy=sgmy*142122.3034;
		 }
		  else
		  {sgmy=0;sgmmy=0;}


		   if(loop3>=0&&loop3<=NNN-1)
		 {  sgmz=sgmmax*(NNN-loop3-1)/NNN*(NNN-loop3-1)/NNN;
		    sgmmz=sgmz*142122.3034;
		 }
		  else if(loop3>=Nz-NNN&&loop3<Nz)
		 {  sgmz=sgmmax*(NNN+loop3-Nz)/NNN*(NNN+loop3-Nz)/NNN;
		    sgmmz=sgmz*142122.3034;
		 }
		  else
		  {sgmz=0;sgmmz=0;}

        if(loop1==Nx/2&&loop2==Ny/2&&loop3>=Nz/2+D/2&&loop3<Nz/2+D/2+DD/2)
		{Ezx[loop1][loop2][loop3]=0;Ezy[loop1][loop2][loop3]=0;}
        else if(loop1==Nx/2&&loop2==Ny/2&&loop3<=Nz/2-D/2&&loop3>Nz/2-D/2-DD/2)
			{Ezx[loop1][loop2][loop3]=0;Ezy[loop1][loop2][loop3]=0;}
		else
		{
          Ezx[loop1][loop2][loop3]=((2*e-sgmx*dt)/(2*e+sgmx*dt))*Ezx[loop1][loop2][loop3]+(2*dt)/(2*e+sgmx*dt)*((Hyz[loop1][loop2][loop3]+Hyx[loop1][loop2][loop3]-Hyz[loop1-1][loop2][loop3]-Hyx[loop1-1][loop2][loop3])/dx);
		  Ezy[loop1][loop2][loop3]=((2*e-sgmy*dt)/(2*e+sgmy*dt))*Ezy[loop1][loop2][loop3]-(2*dt)/(2*e+sgmy*dt)*((Hxy[loop1][loop2][loop3]+Hxz[loop1][loop2][loop3]-Hxy[loop1][loop2-1][loop3]-Hxz[loop1][loop2-1][loop3])/dy);
		}
		   if(loop1>NNN-1&&loop1<Nx-NNN&&loop2>NNN-1&&loop2<Ny-NNN&&loop3>NNN-1&&loop3<Nz-NNN)
1234

⌨️ 快捷键说明

复制代码 Ctrl + C
搜索代码 Ctrl + F
全屏模式 F11
切换主题 Ctrl + Shift + D
显示快捷键 ?
增大字号 Ctrl + =
减小字号 Ctrl + -