h_iter.c

FDTD 模拟 天线 包含 近场－远场 外推程序（时谐场）

/*------------------------------------------
   
     The Iteration Formulae of H Component 

*/

# include "H_Iter.h"

void H_Inner_Iter(MODULE *pModule,double ****Ex,double ****Ey,double ****Ez,
                  double ****Hx,double ****Hy,double ****Hz)
{
  int i,j,k,l,N_k_End;

  int NHxx,NHxy,NHxz;
  int NHyx,NHyy,NHyz;
  int NHzx,NHzy,NHzz;
  double tmp1,tmp2;



  NHxx=Nx+1;  NHxy=Ny;    NHxz=Nz;
  NHyx=Nx;    NHyy=Ny+1;  NHyz=Nz;
  NHzx=Nx;    NHzy=Ny;    NHzz=Nz+1;

  N_k_End=pModule->pLayer[pModule->nLayerNum-1].z_end;

  /* Iteration Formulae of Hx Inner Points */

   for(k=pModule->pLayer[0].z_start; k<N_k_End; k++)
    {
      for(i=1;i<NHxx-1;i++)
       {
         for(j=0;j<NHxy;j++)
          {
            tmp1=(Ey[0][i][j][k+1]-Ey[0][i][j][k])/dhz;
            tmp2=(Ez[0][i][j+1][k]-Ez[0][i][j][k])/dhy;

            Hx[1][i][j][k]=Hx[0][i][j][k]+dt*C*(tmp1-tmp2);

          }
        }
     }      


   /* Iteration Formulae of Hy Inner Points  */

    for(k=pModule->pLayer[0].z_start; k<N_k_End; k++)
    {
      for(i=0;i<NHyx;i++)
       {
         for(j=1;j<NHyy-1;j++)
          {
            tmp1=(Ez[0][i+1][j][k]-Ez[0][i][j][k])/dhx;
            tmp2=(Ex[0][i][j][k+1]-Ex[0][i][j][k])/dhz;

            Hy[1][i][j][k]=Hy[0][i][j][k]+dt*C*(tmp1-tmp2);


          }
        }
     }      

   /*  Iteration Formulae of Hz Inner Points   */

     for(l=0;l<pModule->nLayerNum;l++)
      {
        for(k=pModule->pLayer[l].z_start+1; k<pModule->pLayer[l].z_end;k++)
         {
           for(i=0;i<NHzx;i++)
            {
              for(j=0;j<NHzy;j++)
                {

                  tmp1=(Ex[0][i][j+1][k]-Ex[0][i][j][k])/dhy;
                  tmp2=(Ey[0][i+1][j][k]-Ey[0][i][j][k])/dhx;

                  Hz[1][i][j][k]=Hz[0][i][j][k]+dt*C*(tmp1-tmp2);
               
                 }
             }
         }
     }    


}


/* Judge Hz Point is locate on the Patch

*/


void H_Surface_Iter(MODULE *pModule, double ****Ex, 
                    double ****Ey, double ****Hz)
{

  int i,j,l,m,tmp_k,N_k_End;
  int NHzx,NHzy,NHzz;
  double tmp1,tmp2;
  BOOL Flagc;
  POINT *p;


  NHzx=Nx;    NHzy=Ny;    NHzz=Nz+1;  

  N_k_End=pModule->pLayer[pModule->nLayerNum-1].z_end;


  p=(POINT*) malloc(sizeof(POINT));

  for(l=0;l<pModule->nLayerNum-1;l++)
   {

   tmp_k=pModule->pLayer[l+1].z_start; 

   if(pModule->pLayer[l].nCondPolyNum!=0)  /* Has conductor */  
    {  
     for(i=0;i<NHzx;i++)
     { 
        for(j=0;j<NHzy;j++)
        { 
              p->x=0.5+1.0*i;
              p->y=0.5+1.0*j;
        for(m=0;m<pModule->pLayer[l].nCondPolyNum;m++) /* scan all Conductors */
          {     
               Flagc=OnPolygon(&pModule->pLayer[l].pCondPoly[m],p);   
               if(Flagc==TRUE ) break;
            }
        if(Flagc==TRUE) /* It means that the point p is a conductor point */
           {
              
            Hz[1][i][j][tmp_k]=0.0; 
             
           }
         else
           {
             tmp1=(Ex[0][i][j+1][tmp_k]-Ex[0][i][j][tmp_k])/dhy;
             tmp2=(Ey[0][i+1][j][tmp_k]-Ey[0][i][j][tmp_k])/dhx;

             Hz[1][i][j][tmp_k]=Hz[0][i][j][tmp_k]+dt*C*(tmp1-tmp2);
           }
        }
          
       }
       
      }
    else  /* No Conductors */
     {
       for(i=0;i<NHzx;i++)
         { 
          for(j=0;j<NHzy;j++)
          {
      
             tmp1=(Ex[0][i][j+1][tmp_k]-Ex[0][i][j][tmp_k])/dhy;
             tmp2=(Ey[0][i+1][j][tmp_k]-Ey[0][i][j][tmp_k])/dhx;

             Hz[1][i][j][tmp_k]=Hz[0][i][j][tmp_k]+dt*C*(tmp1-tmp2);
 
           }
          }
       }
   }


    free(p); 

/* -----  Iteration of Hz at the Interface Plane is Over ------- */ 
     
 }


void H_2_Surface_Iter(MODULE *pModule, double ****Ex, 
                    double ****Ey, double ****Hz)
{

  int i,j,tmp_k;
  int NHzx,NHzy,NHzz;
  double tmp1,tmp2;
  BOOL Flagc;
  POINT *p;
  POLYGON *pCondPoly;


  NHzx=Nx;    NHzy=Ny;    NHzz=Nz+1;  


/*  !! Note : Modify the Parameters */ 


  p=(POINT*) malloc(sizeof(POINT));

  pCondPoly=POLYGON_malloc(1);   /*  The Microstrip Parameters  */  

  pCondPoly->nVertex=4;

  pCondPoly->pVertex=POINT_malloc(pCondPoly->nVertex);

  pCondPoly->pVertex[0].x=1.0*pModule->Excite.nFw_start;
  pCondPoly->pVertex[0].y=0.0;

  pCondPoly->pVertex[1].x=1.0*pModule->Excite.nFw_start;
  pCondPoly->pVertex[1].y=1.0*Ny;

  pCondPoly->pVertex[2].x=1.0*pModule->Excite.nFw_end;
  pCondPoly->pVertex[2].y=1.0*Ny;

  pCondPoly->pVertex[3].x=1.0*pModule->Excite.nFw_end;
  pCondPoly->pVertex[3].y=0.0;

   tmp_k=pModule->pLayer[1].z_start; 

     for(i=0;i<NHzx;i++)
     { 
        for(j=0;j<NHzy;j++)
        { 
              p->x=0.5+1.0*i;
              p->y=0.5+1.0*j;
              Flagc=OnPolygon(pCondPoly,p);   

        if(Flagc==TRUE) /* It means that the point p is a conductor point */
           {
              
            Hz[1][i][j][tmp_k]=0.0; 
             
           }
         else
           {
             tmp1=(Ex[0][i][j+1][tmp_k]-Ex[0][i][j][tmp_k])/dhy;
             tmp2=(Ey[0][i+1][j][tmp_k]-Ey[0][i][j][tmp_k])/dhx;

             Hz[1][i][j][tmp_k]=Hz[0][i][j][tmp_k]+dt*C*(tmp1-tmp2);
           }
        }
          
       }

    if (p!=NULL) {
                     free(p); 
                 }  
     if( pCondPoly->pVertex !=NULL)
      {
       free(pCondPoly->pVertex);
      }   

     if( pCondPoly !=NULL)
       { 

          free(pCondPoly);
        }
        

     
 }