"Fossies" - the Fresh Open Source Software Archive

Member "getdp-3.1.0-source/Kernel/GeoNormal.cpp" (29 Dec 2018, 4495 Bytes) of package /linux/privat/getdp-3.1.0-source.tgz:


As a special service "Fossies" has tried to format the requested source page into HTML format using (guessed) C and C++ source code syntax highlighting (style: standard) with prefixed line numbers and code folding option. Alternatively you can here view or download the uninterpreted source code file. For more information about "GeoNormal.cpp" see the Fossies "Dox" file reference documentation and the latest Fossies "Diffs" side-by-side code changes report: 3.0.4_vs_3.1.0.

    1 // GetDP - Copyright (C) 1997-2019 P. Dular and C. Geuzaine, University of Liege
    2 //
    3 // See the LICENSE.txt file for license information. Please report all
    4 // issues on https://gitlab.onelab.info/getdp/getdp/issues.
    5 
    6 #include <string.h>
    7 #include <math.h>
    8 #include "ProData.h"
    9 #include "GeoData.h"
   10 #include "MallocUtils.h"
   11 #include "Message.h"
   12 
   13 #define SQU(a)     ((a)*(a)) 
   14 
   15 extern struct Problem    Problem_S ;
   16 extern struct GeoData  * CurrentGeoData ;
   17 
   18 int  fcmp_NXE(const void * a, const void * b)
   19 {
   20   return 
   21     ((struct Entity2XEntity1 *)a)->Num -
   22     ((struct Entity2XEntity1 *)b)->Num ;
   23 }
   24 
   25 int  fcmp_EXVector(const void * a, const void * b)
   26 {
   27   return 
   28     ((struct EntityXVector *)a)->Num -
   29     ((struct EntityXVector *)b)->Num ;
   30 }
   31 
   32 /*
   33   C'est une maniere un peu naive de creer cette BD. Mais elle a 
   34   l'avantage de permettre une allocation simple (et minimum).
   35 
   36   L'autre possibilite (boucler sur les elemnts) est plus rapide, mais
   37   je ne vois pas bien comment obtenir un cout memoire minimum simplement, 
   38   sans faire de nombreux realloc.
   39 */
   40 
   41 #define MAX_NBR_NXE_INCIDENCE 20 
   42 
   43 static int RegionIndexForNXE = -1;
   44 
   45 void Geo_CreateNodesXElements(int NumNode, int InIndex, 
   46                   int *NbrElements, int **NumElements)
   47 {
   48   struct Entity2XEntity1   NXE, *NXE_P ;
   49   struct Geo_Element      *GeoElement ;
   50   struct Group            *Group_P ;
   51 
   52   int    i, j, tmp[MAX_NBR_NXE_INCIDENCE] ;
   53 
   54   Group_P = (struct Group*)List_Pointer(Problem_S.Group, InIndex);
   55 
   56   if(InIndex != RegionIndexForNXE){
   57     RegionIndexForNXE = InIndex ;
   58     Message::Info("  Generate NodesXElements information for Region '%s'", Group_P->Name);
   59     if(CurrentGeoData->NodesXElements)
   60       Tree_Delete(CurrentGeoData->NodesXElements);
   61     CurrentGeoData->NodesXElements = 
   62       Tree_Create(sizeof(struct Entity2XEntity1), fcmp_NXE) ;
   63   }
   64   
   65   NXE.Num = NumNode ;
   66   
   67   if((NXE_P = (struct Entity2XEntity1*)
   68       Tree_PQuery(CurrentGeoData->NodesXElements, &NXE))) {
   69     *NbrElements = NXE_P->NbrEntities ;
   70     *NumElements = NXE_P->NumEntities ;
   71   }
   72   else{
   73     NXE.NbrEntities = 0 ;
   74     for (i = 0 ; i < Geo_GetNbrGeoElements(); i++) {
   75       GeoElement = Geo_GetGeoElement(i) ;
   76       if (List_Search(Group_P->InitialList, &GeoElement->Region, fcmp_int)){
   77     for(j=0 ; j<GeoElement->NbrNodes ; j++){
   78       if(GeoElement->NumNodes[j] == NumNode){
   79 
   80         /* printf("Adding elm %d to node %d\n", GeoElement->Num, NumNode); */
   81         /* this is to have orientation of elements adjacent to the node 
   82            Only valid for line elemnts !!!!!!! */
   83         
   84         tmp[NXE.NbrEntities] = ((!j)?-1:1) * GeoElement->Num ;
   85         NXE.NbrEntities++ ;
   86       }
   87     }
   88       }
   89     }
   90     NXE.NumEntities = (int*)Malloc(NXE.NbrEntities * sizeof(int)) ;
   91     memcpy(NXE.NumEntities, tmp, NXE.NbrEntities * sizeof(int));
   92     Tree_Add(CurrentGeoData->NodesXElements, &NXE);
   93     *NbrElements = NXE.NbrEntities ;
   94     *NumElements = NXE.NumEntities ;
   95   }
   96 }
   97 
   98 void Geo_CreateNormal(int Type, double *x, double *y, double *z, double *N)
   99 {
  100   double  x1x0, x2x0, y1y0, y2y0, z1z0, z2z0 ;
  101   double  nx, ny, nz, norm ;
  102 
  103   switch (Type) {
  104 
  105   case LINE :
  106   case LINE_2 :
  107     nx = y[1] - y[0] ;
  108     ny = x[0] - x[1] ;
  109     norm = sqrt(SQU(nx)+SQU(ny)) ;      
  110     N[0] = nx / norm ;
  111     N[1] = ny / norm ;
  112     N[2] = 0. ;
  113     break ;
  114 
  115   case TRIANGLE :
  116   case TRIANGLE_2 :
  117   case QUADRANGLE :
  118   case QUADRANGLE_2 :
  119     x1x0 = x[1] - x[0] ;
  120     y1y0 = y[1] - y[0] ;
  121     z1z0 = z[1] - z[0] ;
  122     x2x0 = x[2] - x[0] ; 
  123     y2y0 = y[2] - y[0] ;
  124     z2z0 = z[2] - z[0] ;
  125     nx = y1y0 * z2z0 - z1z0 * y2y0 ;
  126     ny = z1z0 * x2x0 - x1x0 * z2z0 ;
  127     nz = x1x0 * y2y0 - y1y0 * x2x0 ;
  128     norm = sqrt(SQU(nx)+SQU(ny)+SQU(nz)) ;
  129     N[0] = nx/norm ;
  130     N[1] = ny/norm ;
  131     N[2] = nz/norm ;
  132     break ;
  133     
  134   default :
  135     Message::Error("Normal computation not done (yet) for Element Type %d", Type);
  136   }
  137 }
  138 
  139 void Geo_CreateNormalOfElement(struct Geo_Element *GeoElement, double *Normal)
  140 {
  141   struct  EntityXVector  EXV, *EXV_P ;
  142   double  x [NBR_MAX_NODES_IN_ELEMENT] ;
  143   double  y [NBR_MAX_NODES_IN_ELEMENT] ;
  144   double  z [NBR_MAX_NODES_IN_ELEMENT] ;
  145 
  146   EXV.Num = GeoElement->Num ;
  147 
  148   if((EXV_P = (struct EntityXVector*)Tree_PQuery(CurrentGeoData->Normals, &EXV))) {
  149     memcpy(Normal, EXV_P->Vector, 3*sizeof(double));
  150   }
  151   else{
  152     Geo_GetNodesCoordinates(GeoElement->NbrNodes, GeoElement->NumNodes, x, y, z) ;
  153     Geo_CreateNormal(GeoElement->Type, x, y, z, Normal);
  154     memcpy(EXV.Vector, Normal, 3*sizeof(double));
  155     Tree_Add(CurrentGeoData->Normals, &EXV);    
  156   }
  157 }