Salome HOME
Debug on GENERAL_24 and GENERAL_48
[tools/medcoupling.git] / src / INTERP_KERNEL / SplitterTetra.txx
index e36dd0ae55170ba11f1371a016aedc3f76eeeae8..38aab3c233bce06edbf49a475325afa7f9460575 100644 (file)
@@ -1,4 +1,4 @@
-// Copyright (C) 2007-2012  CEA/DEN, EDF R&D
+// Copyright (C) 2007-2013  CEA/DEN, EDF R&D
 //
 // This library is free software; you can redistribute it and/or
 // modify it under the terms of the GNU Lesser General Public
 #include <sstream>
 #include <vector>
 
-/// Smallest volume of the intersecting elements in the transformed space that will be returned as non-zero. 
-/// Since the scale is always the same in the transformed space (the target tetrahedron is unitary), this number is independent of the scale of the meshes.
-#define SPARSE_TRUNCATION_LIMIT 1.0e-14
-
 namespace INTERP_KERNEL
 {
-  
-  /**
-   * Constructor creating object from target cell global number 
-   * 
-   * @param srcMesh     mesh containing the source elements
-   * @param targetMesh  mesh containing the target elements
-   * @param targetCell  global number of the target cell
-   *
-   */
-  /*template<class MyMeshType>
-  SplitterTetra<MyMeshType>::SplitterTetra(const MyMeshType& srcMesh, const MyMeshType& targetMesh, typename MyMeshType::MyConnType targetCell)
-    : _src_mesh(srcMesh), _t(0)
-  {   
-    // get array of corners of target tetraedron
-    const double* tetraCorners[4];
-    for(int i = 0 ; i < 4 ; ++i)
-      tetraCorners[i] = getCoordsOfNode(i, targetCell, targetMesh);
-    // create the affine transform
-    createAffineTransform(tetraCorners);
-    }*/
+  template<class MyMeshType> 
+  const double SplitterTetra<MyMeshType>::SPARSE_TRUNCATION_LIMIT=1.0e-14;
 
   /*!
    * output is expected to be allocated with 24*sizeof(void*) in order to store the 24 tetras.
@@ -81,6 +59,13 @@ namespace INTERP_KERNEL
   }
 
   /**
+   * SplitterTetra class computes for a list of cell ids of a given mesh \a srcMesh (badly named) the intersection with a 
+   * single TETRA4 cell given by \a tetraCorners (of length 4) and \a nodesId (of length 4 too). \a nodedIds is given only to establish
+   * if a partial computation of a triangle has already been performed (to increase performance).
+   *
+   * The \a srcMesh can contain polyhedron cells.
+   * 
+   * 
    * Constructor creating object from the four corners of the tetrahedron.
    *
    * @param srcMesh       mesh containing the source elements
@@ -97,23 +82,9 @@ namespace INTERP_KERNEL
     _coords[6]=tetraCorners[2][0]; _coords[7]=tetraCorners[2][1]; _coords[8]=tetraCorners[2][2];
     _coords[9]=tetraCorners[3][0]; _coords[10]=tetraCorners[3][1]; _coords[11]=tetraCorners[3][2];
     // create the affine transform
-    createAffineTransform(tetraCorners);
+    _t=new TetraAffineTransform(_coords);
   }
-
-  /*!
-   * This contructor is used to build part of 1/24th dual cell of tetraCorners.
-   * @param i is in 0..23 included.
-   * @param nodeId is the id of first node of this in target mesh in C mode.
-   */
-  /*template<class MyMeshType>
-  SplitterTetra<MyMeshType>::SplitterTetra(const MyMeshType& srcMesh, const double** tetraCorners, int i, typename MyMeshType::MyConnType nodeId)
-    : _t(0), _src_mesh(srcMesh), _conn(nodeId)
-  {
-    double *newCoords[4];
-    splitMySelfForDual(tetraCorners,newCoords,i);
-    createAffineTransform(newCoords);
-    }*/
-
+  
   /**
    * Destructor
    *
@@ -566,10 +537,7 @@ namespace INTERP_KERNEL
               counter++;
             }
         }
-      if (counter == 3)
-        return true;
-      else
-        return false;
+      return counter == 3;
   }
 
   /**
@@ -589,7 +557,6 @@ namespace INTERP_KERNEL
    * @param polyNodesNbr number of the nodes of the polygon source face
    * @param polyNodes numbers of the nodes of the polygon source face
    * @param polyCoords coordinates of the nodes of the polygon source face
-   * @param polyCoords coordinates of the nodes of the polygon source face
    * @param dimCaracteristic characteristic size of the meshes containing the triangles
    * @param precision precision for double float data used for comparison
    * @param listOfTetraFacesTreated list of tetra faces treated
@@ -629,9 +596,6 @@ namespace INTERP_KERNEL
         const int globalNodeNum = polyNodes[i];
         if(_nodes.find(globalNodeNum) == _nodes.end())
           {
-            //for(HashMap< int , double* >::iterator iter3=_nodes.begin();iter3!=_nodes.end();iter3++)
-            //  std::cout << (*iter3).first << " ";
-            //std::cout << std::endl << "*** " << globalNodeNum << std::endl;
             calculateNode2(globalNodeNum, polyCoords[i]);
           }
 
@@ -1029,30 +993,7 @@ namespace INTERP_KERNEL
    */
   template<class MyMeshTypeT, class MyMeshTypeS>
   void SplitterTetra2<MyMeshTypeT, MyMeshTypeS>::fiveSplit(const int* const subZone, typename std::vector< SplitterTetra<MyMeshTypeS>* >& tetra)
-  {
-    // Schema according to which the splitting is performed.
-    // Each line represents one tetrahedron. The numbering is as follows :
-    //
-    //          7 ------ 6
-    //         /|       /|
-    //        / |      / |
-    //       3 ------ 2  |
-    //       |  |     |  |
-    //       |  |     |  |
-    //       |  4-----|- 5
-    //       | /      | /
-    //       0 ------ 1
-
-
-    static const int SPLIT_NODES_5[20] = 
-      {
-        0, 1, 5, 2,
-        0, 4, 5, 7,
-        0, 3, 7, 2,
-        5, 6, 7, 2,
-        0, 2, 5, 7
-      };
-    
+  { 
     // create tetrahedra
     for(int i = 0; i < 5; ++i)
       {
@@ -1078,29 +1019,6 @@ namespace INTERP_KERNEL
   template<class MyMeshTypeT, class MyMeshTypeS>
   void SplitterTetra2<MyMeshTypeT, MyMeshTypeS>::sixSplit(const int* const subZone, typename std::vector< SplitterTetra<MyMeshTypeS>* >& tetra)
   {
-    // Schema according to which the splitting is performed.
-    // Each line represents one tetrahedron. The numbering is as follows :
-    //
-    //          7 ------ 6
-    //         /|       /|
-    //        / |      / |
-    //       3 ------ 2  |
-    //       |  |     |  |
-    //       |  |     |  |
-    //       |  4-----|- 5
-    //       | /      | /
-    //       0 ------ 1
-
-    static const int SPLIT_NODES_6[24] = 
-      {
-        0, 1, 5, 6,
-        0, 2, 1, 6,
-        0, 5, 4, 6,
-        0, 4, 7, 6,
-        0, 3, 2, 6,
-        0, 7, 3, 6
-      };
-
     for(int i = 0; i < 6; ++i)
       {
         const double* nodes[4];
@@ -1128,39 +1046,6 @@ namespace INTERP_KERNEL
     // The two nodes of the original mesh cell used in each tetrahedron.
     // The tetrahedra all have nodes (cellCenter, faceCenter, edgeNode1, edgeNode2)
     // For the correspondance of the nodes, see the GENERAL_48_SUB_NODES table in calculateSubNodes
-    static const int TETRA_EDGES[48] = 
-      {
-        // face with center 9
-        0,1,
-        1,5,
-        5,4,
-        4,0,
-        // face with center 10
-        0,1,
-        1,2,
-        2,3,
-        3,0,
-        // face with center 11
-        0,4,
-        4,7,
-        7,3,
-        3,0,
-        // face with center 12
-        1,5,
-        5,6,
-        6,2,
-        2,1,
-        // face with center 13
-        5,6,
-        6,7,
-        7,4,
-        4,5,
-        // face with center 14
-        2,6,
-        6,7,
-        7,3,
-        3,2
-      };
     
     // nodes to use for tetrahedron
     const double* nodes[4];
@@ -1177,8 +1062,8 @@ namespace INTERP_KERNEL
         for(int j = 0; j < 4; ++j)
           {
             const int row = 4*(faceCenterNode - 8) + j;
-            conn[2] = TETRA_EDGES[2*row];
-            conn[3] = TETRA_EDGES[2*row + 1];
+            conn[2] = TETRA_EDGES_GENERAL_24[2*row];
+            conn[3] = TETRA_EDGES_GENERAL_24[2*row + 1];
             nodes[2] = getCoordsOfSubNode(conn[2]);
             nodes[3] = getCoordsOfSubNode(conn[3]);
 
@@ -1199,25 +1084,10 @@ namespace INTERP_KERNEL
    */
   template<class MyMeshTypeT, class MyMeshTypeS>
   void SplitterTetra2<MyMeshTypeT, MyMeshTypeS>::calculateGeneral48Tetra(typename std::vector< SplitterTetra<MyMeshTypeS>* >& tetra)
-  {
-    // Define 8 hexahedral subzones as in Grandy, p449
-    // the values correspond to the nodes that correspond to nodes 1,2,3,4,5,6,7,8 in the subcell
-    // For the correspondance of the nodes, see the GENERAL_48_SUB_NODES table in calculateSubNodes
-    static const int subZones[64] = 
-      {
-        0,8,21,12,9,20,26,22,
-        8,1,13,21,20,10,23,26,
-        12,21,16,3,22,26,25,17,
-        21,13,2,16,26,23,18,25,
-        9,20,26,22,4,11,24,14,
-        20,10,23,26,11,5,15,24,
-        22,26,25,17,14,24,19,7,
-        26,23,18,25,24,15,6,19
-      };
-    
+  { 
     for(int i = 0; i < 8; ++i)
       {
-        sixSplit(&subZones[8*i],tetra);
+        sixSplit(GENERAL_48_SUBZONES+8*i,tetra);
       }
   }
   
@@ -1341,20 +1211,6 @@ namespace INTERP_KERNEL
           {
           case GENERAL_24:
             {
-              // Each sub-node is the barycenter of 4 other nodes.
-              // For the faces, these are on the orignal mesh.
-              // For the barycenter, the four face sub-nodes are used.
-              static const int GENERAL_24_SUB_NODES[28] = 
-                {
-                  0,1,4,5,// sub-node 9 (face)
-                  0,1,2,3,// sub-node 10 (face)
-                  0,3,4,7,// sub-node 11 (face)
-                  1,2,5,6,// sub-node 12 (face)
-                  4,5,6,7,// sub-node 13 (face)
-                  2,3,6,7,// sub-node 14 (face)
-                  8,9,10,11// sub-node 15 (cell)
-                };
-
               for(int i = 0; i < 7; ++i)
                 {
                   double* barycenter = new double[3];
@@ -1366,33 +1222,6 @@ namespace INTERP_KERNEL
 
           case GENERAL_48:
             {
-              // Each sub-node is the barycenter of two other nodes.
-              // For the edges, these lie on the original mesh.
-              // For the faces, these are the edge sub-nodes.
-              // For the cell these are two face sub-nodes.
-              static const int GENERAL_48_SUB_NODES[38] = 
-                {
-                  0,1,   // sub-node 9 (edge)
-                  0,4,   // sub-node 10 (edge)
-                  1,5,   // sub-node 11 (edge)
-                  4,5,   // sub-node 12 (edge)
-                  0,3,   // sub-node 13 (edge)
-                  1,2,   // sub-node 14 (edge)
-                  4,7,   // sub-node 15 (edge)
-                  5,6,   // sub-node 16 (edge)
-                  2,3,   // sub-node 17 (edge)
-                  3,7,   // sub-node 18 (edge)
-                  2,6,   // sub-node 19 (edge)
-                  6,7,   // sub-node 20 (edge)
-                  8,11,  // sub-node 21 (face)
-                  12,13, // sub-node 22 (face)
-                  9,17,  // sub-node 23 (face)
-                  10,18, // sub-node 24 (face)
-                  14,15, // sub-node 25 (face)
-                  16,19, // sub-node 26 (face)
-                  20,25  // sub-node 27 (cell)
-                };
-
               for(int i = 0; i < 19; ++i)
                 {
                   double* barycenter = new double[3];