pos_element.c

cfd求解器使用与gmsh网格的求解

	  PE->u[3] =-1. ; PE->v[3] = 1. ; PE->w[3] = 1. ;
	  POS_CUT_SKIN;
	}
	break ;

      case PRISM :
	PE = Create_PostElement(Index, TRIANGLE, 3, 1) ; /* nodes 1 3 2 */
	PE->NumNodes[0] = GE->NumNodes[0] ;
	PE->NumNodes[1] = GE->NumNodes[2] ;
	PE->NumNodes[2] = GE->NumNodes[1] ;
	PE->u[0] = 0. ; PE->v[0] = 0. ; PE->w[0] =-1. ;
	PE->u[1] = 0. ; PE->v[1] = 1. ; PE->w[1] =-1. ;
	PE->u[2] = 1. ; PE->v[2] = 0. ; PE->w[2] =-1. ;
	POS_CUT_SKIN;
	
	PE = Create_PostElement(Index, TRIANGLE, 3, 1) ; /* nodes 4 5 6 */
	PE->NumNodes[0] = GE->NumNodes[3] ;
	PE->NumNodes[1] = GE->NumNodes[4] ;
	PE->NumNodes[2] = GE->NumNodes[5] ;
	PE->u[0] = 0. ; PE->v[0] = 0. ; PE->w[0] = 1. ;
	PE->u[1] = 1. ; PE->v[1] = 0. ; PE->w[1] = 1. ;
	PE->u[2] = 0. ; PE->v[2] = 1. ; PE->w[2] = 1. ;
	POS_CUT_SKIN;
	
	if(DecomposeInSimplex){
	  PE = Create_PostElement(Index, TRIANGLE, 3, 1) ; /* nodes 1 2 5 */
	  PE->NumNodes[0] = GE->NumNodes[0] ;
	  PE->NumNodes[1] = GE->NumNodes[1] ;
	  PE->NumNodes[2] = GE->NumNodes[4] ;
	  PE->u[0] = 0. ; PE->v[0] = 0. ; PE->w[0] =-1. ;
	  PE->u[1] = 1. ; PE->v[1] = 0. ; PE->w[1] =-1. ;
	  PE->u[2] = 1. ; PE->v[2] = 0. ; PE->w[2] = 1. ;
	  POS_CUT_SKIN;
	
	  PE = Create_PostElement(Index, TRIANGLE, 3, 1) ; /* nodes 1 5 4 */
	  PE->NumNodes[0] = GE->NumNodes[0] ;
	  PE->NumNodes[1] = GE->NumNodes[4] ;
	  PE->NumNodes[2] = GE->NumNodes[3] ;
	  PE->u[0] = 0. ; PE->v[0] = 0. ; PE->w[0] =-1. ;
	  PE->u[1] = 1. ; PE->v[1] = 0. ; PE->w[1] = 1. ;
	  PE->u[2] = 0. ; PE->v[2] = 0. ; PE->w[2] = 1. ;
	  POS_CUT_SKIN;
	  
	  PE = Create_PostElement(Index, TRIANGLE, 3, 1) ; /* nodes 1 6 3 */
	  PE->NumNodes[0] = GE->NumNodes[0] ;
	  PE->NumNodes[1] = GE->NumNodes[5] ;
	  PE->NumNodes[2] = GE->NumNodes[2] ;
	  PE->u[0] = 0. ; PE->v[0] = 0. ; PE->w[0] =-1. ;
	  PE->u[1] = 0. ; PE->v[1] = 1. ; PE->w[1] = 1. ;
	  PE->u[2] = 0. ; PE->v[2] = 1. ; PE->w[2] =-1. ;
	  POS_CUT_SKIN;
	  
	  PE = Create_PostElement(Index, TRIANGLE, 3, 1) ; /* nodes 1 4 6 */
	  PE->NumNodes[0] = GE->NumNodes[0] ;
	  PE->NumNodes[1] = GE->NumNodes[3] ;
	  PE->NumNodes[2] = GE->NumNodes[5] ;
	  PE->u[0] = 0. ; PE->v[0] = 0. ; PE->w[0] =-1. ;
	  PE->u[1] = 0. ; PE->v[1] = 0. ; PE->w[1] = 1. ;
	  PE->u[2] = 0. ; PE->v[2] = 1. ; PE->w[2] = 1. ;
	  POS_CUT_SKIN;
	
	  PE = Create_PostElement(Index, TRIANGLE, 3, 1) ; /* nodes 2 3 5 */
	  PE->NumNodes[0] = GE->NumNodes[1] ;
	  PE->NumNodes[1] = GE->NumNodes[2] ;
	  PE->NumNodes[2] = GE->NumNodes[4] ;
	  PE->u[0] = 1. ; PE->v[0] = 0. ; PE->w[0] =-1. ;
	  PE->u[1] = 0. ; PE->v[1] = 1. ; PE->w[1] =-1. ;
	  PE->u[2] = 1. ; PE->v[2] = 0. ; PE->w[2] = 1. ;
	  POS_CUT_SKIN;
	
	  PE = Create_PostElement(Index, TRIANGLE, 3, 1) ; /* nodes 3 6 5 */
	  PE->NumNodes[0] = GE->NumNodes[2] ;
	  PE->NumNodes[1] = GE->NumNodes[5] ;
	  PE->NumNodes[2] = GE->NumNodes[4] ;
	  PE->u[0] = 0. ; PE->v[0] = 1. ; PE->w[0] =-1. ;
	  PE->u[1] = 0. ; PE->v[1] = 1. ; PE->w[1] = 1. ;      
	  PE->u[2] = 1. ; PE->v[2] = 0. ; PE->w[2] = 1. ;
	  POS_CUT_SKIN;
	}
	else{
	  PE = Create_PostElement(Index, QUADRANGLE, 4, 1) ; /* nodes 1 2 5 4 */
	  PE->NumNodes[0] = GE->NumNodes[0] ;
	  PE->NumNodes[1] = GE->NumNodes[1] ;
	  PE->NumNodes[2] = GE->NumNodes[4] ;
	  PE->NumNodes[3] = GE->NumNodes[3] ;
	  PE->u[0] = 0. ; PE->v[0] = 0. ; PE->w[0] =-1. ;
	  PE->u[1] = 1. ; PE->v[1] = 0. ; PE->w[1] =-1. ;
	  PE->u[2] = 1. ; PE->v[2] = 0. ; PE->w[2] = 1. ;
	  PE->u[3] = 0. ; PE->v[3] = 0. ; PE->w[3] = 1. ;
	  POS_CUT_SKIN;
	  
	  PE = Create_PostElement(Index, QUADRANGLE, 4, 1) ; /* nodes 1 4 6 3 */
	  PE->NumNodes[0] = GE->NumNodes[0] ;
	  PE->NumNodes[1] = GE->NumNodes[3] ;
	  PE->NumNodes[2] = GE->NumNodes[5] ;
	  PE->NumNodes[3] = GE->NumNodes[2] ;
	  PE->u[0] = 0. ; PE->v[0] = 0. ; PE->w[0] =-1. ;
	  PE->u[1] = 0. ; PE->v[1] = 0. ; PE->w[1] = 1. ;
	  PE->u[2] = 0. ; PE->v[2] = 1. ; PE->w[2] = 1. ;
	  PE->u[3] = 0. ; PE->v[3] = 1. ; PE->w[3] =-1. ;
	  POS_CUT_SKIN;
	  
	  PE = Create_PostElement(Index, QUADRANGLE, 4, 1) ; /* nodes 2 3 6 5 */
	  PE->NumNodes[0] = GE->NumNodes[1] ;
	  PE->NumNodes[1] = GE->NumNodes[2] ;
	  PE->NumNodes[2] = GE->NumNodes[5] ;
	  PE->NumNodes[3] = GE->NumNodes[4] ;
	  PE->u[0] = 1. ; PE->v[0] = 0. ; PE->w[0] =-1. ;
	  PE->u[1] = 0. ; PE->v[1] = 1. ; PE->w[1] =-1. ;
	  PE->u[2] = 0. ; PE->v[2] = 1. ; PE->w[2] = 1. ;
	  PE->u[3] = 1. ; PE->v[3] = 0. ; PE->w[3] = 1. ;
	  POS_CUT_SKIN;
	}
	break ;
	
      case PYRAMID :
	PE = Create_PostElement(Index, TRIANGLE, 3, 1) ; /* nodes 1 2 5 */
	PE->NumNodes[0] = GE->NumNodes[0] ;
	PE->NumNodes[1] = GE->NumNodes[1] ;
	PE->NumNodes[2] = GE->NumNodes[4] ;
#if defined(NEW_PYRAMIDS)
	PE->u[0] = 0. ; PE->v[0] = 0. ; PE->w[0] = 0. ;
	PE->u[1] = 1. ; PE->v[1] = 0. ; PE->w[1] = 0. ;
	PE->u[2] = 0. ; PE->v[2] = 0. ; PE->w[2] = 1. ;
#else
	PE->u[0] =-1. ; PE->v[0] =-1. ; PE->w[0] = 0. ;
	PE->u[1] = 1. ; PE->v[1] =-1. ; PE->w[1] = 0. ;
	PE->u[2] = 0. ; PE->v[2] = 0. ; PE->w[2] = 1. ;
#endif
	POS_CUT_SKIN;

	PE = Create_PostElement(Index, TRIANGLE, 3, 1) ; /* nodes 2 3 5 */
	PE->NumNodes[0] = GE->NumNodes[1] ;
	PE->NumNodes[1] = GE->NumNodes[2] ;
	PE->NumNodes[2] = GE->NumNodes[4] ;
#if defined(NEW_PYRAMIDS)
	PE->u[0] = 1. ; PE->v[0] = 0. ; PE->w[0] = 0. ;
	PE->u[1] = 1. ; PE->v[1] = 1. ; PE->w[1] = 0. ;
	PE->u[2] = 0. ; PE->v[2] = 0. ; PE->w[2] = 1. ;
#else
	PE->u[0] = 1. ; PE->v[0] =-1. ; PE->w[0] = 0. ;
	PE->u[1] = 1. ; PE->v[1] = 1. ; PE->w[1] = 0. ;
	PE->u[2] = 0. ; PE->v[2] = 0. ; PE->w[2] = 1. ;
#endif
	POS_CUT_SKIN;

	PE = Create_PostElement(Index, TRIANGLE, 3, 1) ; /* nodes 3 4 5 */
	PE->NumNodes[0] = GE->NumNodes[2] ;
	PE->NumNodes[1] = GE->NumNodes[3] ;
	PE->NumNodes[2] = GE->NumNodes[4] ;
#if defined(NEW_PYRAMIDS)
	PE->u[0] = 1. ; PE->v[0] = 1. ; PE->w[0] = 0. ;
	PE->u[1] = 0. ; PE->v[1] = 1. ; PE->w[1] = 0. ;
	PE->u[2] = 0. ; PE->v[2] = 0. ; PE->w[2] = 1. ;
#else
	PE->u[0] = 1. ; PE->v[0] = 1. ; PE->w[0] = 0. ;
	PE->u[1] =-1. ; PE->v[1] = 1. ; PE->w[1] = 0. ;
	PE->u[2] = 0. ; PE->v[2] = 0. ; PE->w[2] = 1. ;
#endif
	POS_CUT_SKIN;

	PE = Create_PostElement(Index, TRIANGLE, 3, 1) ; /* nodes 4 1 5 */
	PE->NumNodes[0] = GE->NumNodes[3] ;
	PE->NumNodes[1] = GE->NumNodes[0] ;
	PE->NumNodes[2] = GE->NumNodes[4] ;
#if defined(NEW_PYRAMIDS)
	PE->u[0] = 0. ; PE->v[0] = 1. ; PE->w[0] = 0. ;
	PE->u[1] = 0. ; PE->v[1] = 0. ; PE->w[1] = 0. ;
	PE->u[2] = 0. ; PE->v[2] = 0. ; PE->w[2] = 1. ;
#else
	PE->u[0] =-1. ; PE->v[0] = 1. ; PE->w[0] = 0. ;
	PE->u[1] =-1. ; PE->v[1] =-1. ; PE->w[1] = 0. ;
	PE->u[2] = 0. ; PE->v[2] = 0. ; PE->w[2] = 1. ;
#endif
	POS_CUT_SKIN;
	
	if(DecomposeInSimplex){
	  PE = Create_PostElement(Index, TRIANGLE, 3, 1) ; /* nodes 1 3 2 */
	  PE->NumNodes[0] = GE->NumNodes[0] ;
	  PE->NumNodes[1] = GE->NumNodes[2] ;
	  PE->NumNodes[2] = GE->NumNodes[1] ;
#if defined(NEW_PYRAMIDS)
	  PE->u[0] = 0. ; PE->v[0] = 0. ; PE->w[0] = 0. ;
	  PE->u[1] = 1. ; PE->v[1] = 1. ; PE->w[1] = 0. ;
	  PE->u[2] = 1. ; PE->v[2] = 0. ; PE->w[2] = 0. ;
#else
	  PE->u[0] =-1. ; PE->v[0] =-1. ; PE->w[0] = 0. ;
	  PE->u[1] = 1. ; PE->v[1] = 1. ; PE->w[1] = 0. ;
	  PE->u[2] = 1. ; PE->v[2] =-1. ; PE->w[2] = 0. ;
#endif
	  POS_CUT_SKIN;
	  
	  PE = Create_PostElement(Index, TRIANGLE, 3, 1) ; /* nodes 1 4 3 */
	  PE->NumNodes[0] = GE->NumNodes[0] ;
	  PE->NumNodes[1] = GE->NumNodes[3] ;
	  PE->NumNodes[2] = GE->NumNodes[2] ;
#if defined(NEW_PYRAMIDS)
	  PE->u[0] = 0. ; PE->v[0] = 0. ; PE->w[0] = 0. ;
	  PE->u[1] = 0. ; PE->v[1] = 1. ; PE->w[1] = 0. ;
	  PE->u[2] = 1. ; PE->v[2] = 1. ; PE->w[2] = 0. ;
#else
	  PE->u[0] =-1. ; PE->v[0] =-1. ; PE->w[0] = 0. ;
	  PE->u[1] =-1. ; PE->v[1] = 1. ; PE->w[1] = 0. ;
	  PE->u[2] = 1. ; PE->v[2] = 1. ; PE->w[2] = 0. ;
#endif
	  POS_CUT_SKIN;	
	}
	else{
	  PE = Create_PostElement(Index, QUADRANGLE, 4, 1) ; /* nodes 1 4 3 2 */
	  PE->NumNodes[0] = GE->NumNodes[0] ;
	  PE->NumNodes[1] = GE->NumNodes[3] ;
	  PE->NumNodes[2] = GE->NumNodes[2] ;
	  PE->NumNodes[3] = GE->NumNodes[1] ;
#if defined(NEW_PYRAMIDS)
	  PE->u[0] = 0. ; PE->v[0] = 0. ; PE->w[0] = 0. ;
	  PE->u[1] = 0. ; PE->v[1] = 1. ; PE->w[1] = 0. ;
	  PE->u[2] = 1. ; PE->v[2] = 1. ; PE->w[2] = 0. ;
	  PE->u[3] = 1. ; PE->v[3] = 0. ; PE->w[3] = 0. ;
#else
	  PE->u[0] =-1. ; PE->v[0] =-1. ; PE->w[0] = 0. ;
	  PE->u[1] =-1. ; PE->v[1] = 1. ; PE->w[1] = 0. ;
	  PE->u[2] = 1. ; PE->v[2] = 1. ; PE->w[2] = 0. ;
	  PE->u[3] = 1. ; PE->v[3] =-1. ; PE->w[3] = 0. ;
#endif
	  POS_CUT_SKIN;
	}
	break ;

      }

    }

  }

  GetDP_End ;
}

#undef POS_CUT_FILL
#undef POS_CUT_SKIN



/* ------------------------------------------------------------------------ */
/*  S o r t B y C o n n e c t i v i t y                                     */
/* ------------------------------------------------------------------------ */

int Compare_PostElement_Node(struct PostElement * PE1, int n1,
			     struct PostElement * PE2, int n2){
  double TOL=Current.GeoData->CharacteristicLength * 1.e-8;
  if ( (fabs(PE1->x[n1] - PE2->x[n2]) < TOL) &&
       (fabs(PE1->y[n1] - PE2->y[n2]) < TOL) &&
       (fabs(PE1->z[n1] - PE2->z[n2]) < TOL) )
    return 1 ;
  return 0 ;
}

void Sort_PostElement_Connectivity(List_T *PostElement_L){
  int ii, jj, start, end, iPost, NbrPost ;
  struct PostElement *PE, *PE2 ;

  GetDP_Begin("Sort_PostElement_Connectivity");

  NbrPost = List_Nbr(PostElement_L) ;

  /* 
     u[0] = 1 if the element is in the ordered list, with natural orientation
           -1 if the element is in the ordered list, but with opposite orientation 
            0 if the element is not in the list
     v[0] = relative index (to the first element) in the ordered list
  */

  for(ii = 0 ; ii < NbrPost ; ii++){
    PE = *(struct PostElement**)List_Pointer(PostElement_L, ii);
    if(PE->NbrNodes != 2)
      Msg(GERROR, "Connectivity sorting impossible for %d-noded elements",
	  PE->NbrNodes) ;
    PE->u[0] = 0. ;
  }

  PE = *(struct PostElement**)List_Pointer(PostElement_L, 0);
  PE->u[0] = 1. ; PE->v[0] = 0. ;

  iPost = 1 ;
  while(iPost < NbrPost){
    for(ii = 0 ; ii < NbrPost ; ii++){
      PE = *(struct PostElement**)List_Pointer(PostElement_L, ii);
      if(PE->u[0]){
	if(PE->u[0] > 0){
	  start = 0 ; end = 1 ;
	}
	else{
	  start = 1 ; end = 0 ;
	}
	for(jj = 0 ; jj < NbrPost ; jj++){
	  if(jj != ii){
	    PE2 = *(struct PostElement**)List_Pointer(PostElement_L, jj);
	    if(!PE2->u[0]){
	      if(Compare_PostElement_Node(PE, end, PE2, 0)){
		PE2->u[0] = 1. ; PE2->v[0] = PE->v[0] + 1. ; iPost++ ;
	      }
	      else if (Compare_PostElement_Node(PE, start, PE2, 0)){
		PE2->u[0] = -1. ; PE2->v[0] = PE->v[0] - 1.  ; iPost++ ;
	      }
	      else if (Compare_PostElement_Node(PE, start, PE2, 1)){
		PE2->u[0] = 1. ; PE2->v[0] = PE->v[0] - 1. ; iPost++ ;
	      }
	      else if (Compare_PostElement_Node(PE, end, PE2, 1)){
		PE2->u[0] = -1. ; PE2->v[0] = PE->v[0] + 1. ; iPost++ ;
	      }
	    }
	  }
	}
      }
    }
    List_Sort(PostElement_L, fcmp_PostElement_absu0) ;
  }
  List_Sort(PostElement_L, fcmp_PostElement_v0) ;

  GetDP_End ;
}