Salome HOME
MEDCoupling1SGTUMesh::sortHexa8EachOther
[tools/medcoupling.git] / src / INTERP_KERNEL / SplitterTetra.hxx
index a5bec7efb4c585e4162c72c0c4f6e0ac9489eecb..a8be6303ac33ffbe62976b8ce3f92ad6d9574414 100644 (file)
-//  Copyright (C) 2007-2008  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
-//  License as published by the Free Software Foundation; either
-//  version 2.1 of the License.
+// This library is free software; you can redistribute it and/or
+// modify it under the terms of the GNU Lesser General Public
+// License as published by the Free Software Foundation; either
+// version 2.1 of the License.
 //
-//  This library is distributed in the hope that it will be useful,
-//  but WITHOUT ANY WARRANTY; without even the implied warranty of
-//  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
-//  Lesser General Public License for more details.
+// This library is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+// Lesser General Public License for more details.
 //
-//  You should have received a copy of the GNU Lesser General Public
-//  License along with this library; if not, write to the Free Software
-//  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
+// You should have received a copy of the GNU Lesser General Public
+// License along with this library; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
 //
-//  See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
 //
+
 #ifndef __SPLITTERTETRA_HXX__
 #define __SPLITTERTETRA_HXX__
 
+#include "INTERPKERNELDefines.hxx"
 #include "TransformedTriangle.hxx"
 #include "TetraAffineTransform.hxx"
+#include "InterpolationOptions.hxx"
+#include "InterpKernelException.hxx"
+#include "InterpKernelHashMap.hxx"
+#include "VectorUtils.hxx"
 
-#include <assert.h>
-#include <vector>
 #include <functional>
+#include <vector>
+#include <cassert>
 #include <map>
-#ifdef WIN32
-# include <hash_map>
-#else
-# include <ext/hash_map>
-#endif
-
-#ifndef WIN32
-using __gnu_cxx::hash_map;
-#else
-using stdext::hash_map;
-using stdext::hash_compare;
-#endif
+#include <set>
 
 namespace INTERP_KERNEL
 {
+  // 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] = /* WHY not all well oriented ???? */
+    {
+      0, 1, 5, 2,
+      0, 4, 5, 7,
+      0, 3, 7, 2,
+      5, 6, 7, 2,
+      0, 2, 5, 7
+    };
+
+  static const int SPLIT_NODES_5_WO[20] = /* WO for well oriented !!! normals of 3 first points are OUTSIDE the TETRA4 */
+    {
+      0, 5, 1, 2,
+      0, 4, 5, 7,
+      0, 3, 7, 2,
+      5, 7, 6, 2,
+      0, 5, 2, 7
+    };
+
+  static const int SPLIT_NODES_6[24] = /* WHY all badly oriented ???? */
+    {
+      0, 1, 5, 6,
+      0, 2, 1, 6,
+      0, 5, 4, 6,
+      0, 4, 7, 6,
+      0, 3, 2, 6,
+      0, 7, 3, 6
+    };
+  
+  static const int SPLIT_NODES_6_WO[24] = /* WO for well oriented !!! normals of 3 first points are OUTSIDE the TETRA4 */
+    {
+      0, 5, 1, 6,
+      0, 1, 2, 6,
+      0, 4, 5, 6,
+      0, 7, 4, 6,
+      0, 2, 3, 6,
+      0, 3, 7, 6
+    };
+  
+  // 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 8  (face)
+      0,1,2,3,// sub-node 9  (face)
+      0,3,4,7,// sub-node 10 (face)
+      1,2,5,6,// sub-node 11 (face)
+      4,5,6,7,// sub-node 12 (face)
+      2,3,6,7,// sub-node 13 (face)
+      8,9,10,11// sub-node 14 (cell)
+    };
+
+  static const int GENERAL_24_SUB_NODES_WO[28] = 
+    {
+      0,4,5,1,// sub-node 8  (face)
+      0,1,2,3,// sub-node 9  (face)
+      0,3,7,4,// sub-node 10 (face)
+      1,5,6,2,// sub-node 11 (face)
+      4,7,6,5,// sub-node 12 (face)
+      2,6,7,3,// sub-node 13 (face)
+      8,9,10,11// sub-node 14 (cell)
+    };
+  
+  static const int TETRA_EDGES_GENERAL_24[48] = 
+    {
+      // face with center 8
+      0,1,
+      1,5,
+      5,4,
+      4,0,
+      // face with center 9
+      0,1,
+      1,2,
+      2,3,
+      3,0,
+      // face with center 10
+      0,4,
+      4,7,
+      7,3,
+      3,0,
+      // face with center 11
+      1,5,
+      5,6,
+      6,2,
+      2,1,
+      // face with center 12
+      5,6,
+      6,7,
+      7,4,
+      4,5,
+      // face with center 13
+      2,6,
+      6,7,
+      7,3,
+      3,2
+    };
+  
+  // 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 8 (edge)
+      0,4,   // sub-node 9 (edge)
+      1,5,   // sub-node 10 (edge)
+      4,5,   // sub-node 11 (edge)
+      0,3,   // sub-node 12 (edge)
+      1,2,   // sub-node 13 (edge)
+      4,7,   // sub-node 14 (edge)
+      5,6,   // sub-node 15 (edge)
+      2,3,   // sub-node 16 (edge)
+      3,7,   // sub-node 17 (edge)
+      2,6,   // sub-node 18 (edge)
+      6,7,   // sub-node 19 (edge)
+      8,11,  // sub-node 20 (face)
+      12,13, // sub-node 21 (face)
+      9,17,  // sub-node 22 (face)
+      10,18, // sub-node 23 (face)
+      14,15, // sub-node 24 (face)
+      16,19, // sub-node 25 (face)
+      20,25  // sub-node 26 (cell)
+    };
+
+  // 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 GENERAL_48_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
+    };
+
+  static const int GENERAL_48_SUBZONES_2[64] = 
+    {
+      0,-1,-14,-5,-2,-13,-19,-15,
+      -1,1,-6,-14,-13,-3,-16,-19,
+      -5,-14,-9,3,-15,-19,-18,-10,
+      -14,-6,2,-9,-19,-16,-11,-18,
+      -2,-13,-19,-15,4,-4,-17,-7,
+      -13,-3,-16,-19,-4,5,-8,-17,
+      -15,-19,-18,-10,-7,-17,-12,7,
+      -19,-16,-11,-18,-17,-8,6,-12};
+
+  void SplitHexa8IntoTetras(SplittingPolicy policy, const int *nodalConnBg, const int *nodalConnEnd, const double *coords,
+                            std::vector<int>& tetrasNodalConn, std::vector<double>& addCoords);
+  
+  INTERPKERNEL_EXPORT void SplitIntoTetras(SplittingPolicy policy, NormalizedCellType gt, const int *nodalConnBg, const int *nodalConnEnd, const double *coords,
+                                           std::vector<int>& tetrasNodalConn, std::vector<double>& addCoords);
+  
   /**
    * \brief Class representing a triangular face, used as key in caching hash map in SplitterTetra.
    *
@@ -60,7 +224,7 @@ namespace INTERP_KERNEL
      */
     TriangleFaceKey(int node1, int node2, int node3)
     {
-      sort3Ints(_nodes, node1, node2, node3);
+      Sort3Ints(_nodes, node1, node2, node3);
       _hashVal = ( _nodes[0] + _nodes[1] + _nodes[2] ) % 29;
     }
 
@@ -76,6 +240,28 @@ namespace INTERP_KERNEL
       return _nodes[0] == key._nodes[0] && _nodes[1] == key._nodes[1] && _nodes[2] == key._nodes[2];
     }
 
+    /**
+     * Less than operator.
+     *
+     * @param   key  TriangleFaceKey with which to compare
+     * @return  true if this object has the three nodes less than the nodes of the key object, false if not
+     */
+    bool operator<(const TriangleFaceKey& key) const
+    {
+      for (int i = 0; i < 3; ++i)
+        {
+          if (_nodes[i] < key._nodes[i])
+            {
+              return true;
+            }
+          else if (_nodes[i] > key._nodes[i])
+            {
+              return false;
+            }
+        }
+      return false;
+    }
+
     /**
      * Returns a hash value for the object, based on its three nodes.
      * This value is not unique for each face.
@@ -86,15 +272,8 @@ namespace INTERP_KERNEL
     {
       return _hashVal;
     }
-
-#ifdef WIN32
-    operator size_t () const
-    {
-      return _hashVal;
-    }
-#endif
      
-    inline void sort3Ints(int* sorted, int node1, int node2, int node3);
+    inline static void Sort3Ints(int* sorted, int node1, int node2, int node3);
 
   private:
     /// global numbers of the three nodes, sorted in ascending order
@@ -112,7 +291,7 @@ namespace INTERP_KERNEL
    * @param x2   second integer
    * @param x3   third integer
    */
-  inline void TriangleFaceKey::sort3Ints(int* sorted, int x1, int x2, int x3)
+  inline void TriangleFaceKey::Sort3Ints(int* sorted, int x1, int x2, int x3)
   {
     if(x1 < x2)
       {
@@ -149,20 +328,13 @@ namespace INTERP_KERNEL
           }
       }
   }
-}
-#ifndef WIN32
-namespace __gnu_cxx
-{
-
 
   /**
-   * \brief Template specialization of __gnu_cxx::hash<T> function object for use with a __gnu_cxx::hash_map 
+   * \brief Template specialization of INTERP_KERNEL::hash<T> function object for use with a 
    * with TriangleFaceKey as key class.
    * 
    */
-  template<>
-  class hash<INTERP_KERNEL::TriangleFaceKey>
-
+  template<> class hash<INTERP_KERNEL::TriangleFaceKey>
   {
   public:
     /**
@@ -177,20 +349,9 @@ namespace __gnu_cxx
     }
   };
 }
-#else
-  struct TriangleFaceKeyComparator
-  {
-    bool operator()(const INTERP_KERNEL::TriangleFaceKey& key1,
-                    const INTERP_KERNEL::TriangleFaceKey& key2 ) const
-    {
-      return key1.hashVal() < key2.hashVal();
-    }
-  };
-#endif
 
 namespace INTERP_KERNEL
 {
-
   /** 
    * \brief Class calculating the volume of intersection between a tetrahedral target element and
    * source elements with triangular or quadratilateral faces.
@@ -203,9 +364,21 @@ namespace INTERP_KERNEL
     
     SplitterTetra(const MyMeshType& srcMesh, const double** tetraCorners, const typename MyMeshType::MyConnType *nodesId);
 
+    SplitterTetra(const MyMeshType& srcMesh, const double tetraCorners[12], const int *conn = 0);
+
     ~SplitterTetra();
 
-    double intersectSourceCell(typename MyMeshType::MyConnType srcCell);
+    double intersectSourceCell(typename MyMeshType::MyConnType srcCell, double* baryCentre=0);
+    double intersectSourceFace(const NormalizedCellType polyType,
+                               const int polyNodesNbr,
+                               const int *const polyNodes,
+                               const double *const *const polyCoords,
+                               const double dimCaracteristic,
+                               const double precision,
+                               std::multiset<TriangleFaceKey>& listOfTetraFacesTreated,
+                               std::set<TriangleFaceKey>& listOfTetraFacesColinear);
+
+    double intersectTetra(const double** tetraCorners);
 
     typename MyMeshType::MyConnType getId(int id) { return _conn[id]; }
     
@@ -213,13 +386,23 @@ namespace INTERP_KERNEL
 
     void splitMySelfForDual(double* output, int i, typename MyMeshType::MyConnType& nodeId);
 
+    void clearVolumesCache();
+
   private:
-    // member functions
-    inline void createAffineTransform(const double** corners);
-    inline void checkIsOutside(const double* pt, bool* isOutside) const;
+    inline static void CheckIsOutside(const double* pt, bool* isOutside, const double errTol = DEFAULT_ABS_TOL);
+    inline static void CheckIsStrictlyOutside(const double* pt, bool* isStrictlyOutside, const double errTol = DEFAULT_ABS_TOL);
     inline void calculateNode(typename MyMeshType::MyConnType globalNodeNum);
+    inline void calculateNode2(typename MyMeshType::MyConnType globalNodeNum, const double* node);
     inline void calculateVolume(TransformedTriangle& tri, const TriangleFaceKey& key);
-        
+    inline void calculateSurface(TransformedTriangle& tri, const TriangleFaceKey& key);
+
+    static inline bool IsFacesCoplanar(const double *const planeNormal, const double planeConstant,
+                                const double *const *const coordsFace, const double precision);
+    static inline double CalculateIntersectionSurfaceOfCoplanarTriangles(const double *const planeNormal,
+                                                                         const double planeConstant,
+                                                                         const double *const p1, const double *const p2, const double *const p3,
+                                                                         const double *const p4, const double *const p5, const double *const p6,
+                                                                         const double dimCaracteristic, const double precision);
 
     /// disallow copying
     SplitterTetra(const SplitterTetra& t);
@@ -231,15 +414,11 @@ namespace INTERP_KERNEL
     /// affine transform associated with this target element
     TetraAffineTransform* _t;
     
-    /// hash_map relating node numbers to transformed nodes, used for caching
-    hash_map< int , double* > _nodes;
+    /// HashMap relating node numbers to transformed nodes, used for caching
+    HashMap< int , double* > _nodes;
     
-    /// hash_map relating triangular faces to calculated volume contributions, used for caching
-    hash_map< TriangleFaceKey, double
-#ifdef WIN32
-        , hash_compare<TriangleFaceKey,TriangleFaceKeyComparator> 
-#endif
-    > _volumes;
+    /// HashMap relating triangular faces to calculated volume contributions, used for caching
+    HashMap< TriangleFaceKey, double > _volumes;
 
     /// reference to the source mesh
     const MyMeshType& _src_mesh;
@@ -248,21 +427,12 @@ namespace INTERP_KERNEL
     typename MyMeshType::MyConnType _conn[4];
 
     double _coords[12];
+    
+    /// 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.
+    static const double SPARSE_TRUNCATION_LIMIT;
   };
 
-  /**
-   * Creates the affine transform _t from the corners of the tetrahedron. Used by the constructors
-   *
-   * @param corners  double*[4] array containing pointers to four double[3] arrays with the 
-   *                 coordinates of the corners of the tetrahedron
-   */
-  template<class MyMeshType>
-  inline void SplitterTetra<MyMeshType>::createAffineTransform(const double** corners)
-  {
-    // create AffineTransform from tetrahedron
-    _t = new TetraAffineTransform( corners );
-  }
-
   /**
    * Function used to filter out elements by checking if they belong to one of the halfspaces
    * x <= 0, x >= 1, y <= 0, y >= 1, z <= 0, z >= 1, (indexed 0 - 7). The function updates an array of boolean variables
@@ -273,18 +443,31 @@ namespace INTERP_KERNEL
    * @param isOutside bool[8] which indicate the results of earlier checks. 
    */
   template<class MyMeshType>
-  inline void SplitterTetra<MyMeshType>::checkIsOutside(const double* pt, bool* isOutside) const
+  inline void SplitterTetra<MyMeshType>::CheckIsOutside(const double* pt, bool* isOutside, const double errTol)
   {
-    isOutside[0] = isOutside[0] && (pt[0] <= 0.0);
-    isOutside[1] = isOutside[1] && (pt[0] >= 1.0);
-    isOutside[2] = isOutside[2] && (pt[1] <= 0.0);
-    isOutside[3] = isOutside[3] && (pt[1] >= 1.0);
-    isOutside[4] = isOutside[4] && (pt[2] <= 0.0);
-    isOutside[5] = isOutside[5] && (pt[2] >= 1.0);
-    isOutside[6] = isOutside[6] && (1.0 - pt[0] - pt[1] - pt[2] <= 0.0);
-    isOutside[7] = isOutside[7] && (1.0 - pt[0] - pt[1] - pt[2] >= 1.0);
+    isOutside[0] = isOutside[0] && (pt[0] < errTol);
+    isOutside[1] = isOutside[1] && (pt[0] > (1.0-errTol) );
+    isOutside[2] = isOutside[2] && (pt[1] < errTol);
+    isOutside[3] = isOutside[3] && (pt[1] > (1.0-errTol));
+    isOutside[4] = isOutside[4] && (pt[2] < errTol);
+    isOutside[5] = isOutside[5] && (pt[2] > (1.0-errTol));
+    isOutside[6] = isOutside[6] && (1.0 - pt[0] - pt[1] - pt[2] < errTol);
+    isOutside[7] = isOutside[7] && (1.0 - pt[0] - pt[1] - pt[2] > (1.0-errTol) );
   }
   
+  template<class MyMeshType>
+  inline void SplitterTetra<MyMeshType>::CheckIsStrictlyOutside(const double* pt, bool* isStrictlyOutside, const double errTol)
+  {
+    isStrictlyOutside[0] = isStrictlyOutside[0] && (pt[0] < -errTol);
+    isStrictlyOutside[1] = isStrictlyOutside[1] && (pt[0] > (1.0 + errTol));
+    isStrictlyOutside[2] = isStrictlyOutside[2] && (pt[1] < -errTol);
+    isStrictlyOutside[3] = isStrictlyOutside[3] && (pt[1] > (1.0 + errTol));
+    isStrictlyOutside[4] = isStrictlyOutside[4] && (pt[2] < -errTol);
+    isStrictlyOutside[5] = isStrictlyOutside[5] && (pt[2] > (1.0 + errTol));
+    isStrictlyOutside[6] = isStrictlyOutside[6] && (1.0 - pt[0] - pt[1] - pt[2] < -errTol);
+    isStrictlyOutside[7] = isStrictlyOutside[7] && (1.0 - pt[0] - pt[1] - pt[2] > (1.0 + errTol));
+  }
+
   /**
    * Calculates the transformed node with a given global node number.
    * Gets the coordinates for the node in _src_mesh with the given global number and applies TetraAffineTransform
@@ -304,6 +487,26 @@ namespace INTERP_KERNEL
     _nodes[globalNodeNum] = transformedNode;
   }
 
+
+  /**
+   * Calculates the transformed node with a given global node number.
+   * Applies TetraAffineTransform * _t to it.
+   * Stores the result in the cache _nodes. The non-existence of the node in _nodes should be verified before * calling.
+   * The only difference with the previous method calculateNode is that the coordinates of the node are passed in arguments
+   * and are not recalculated in order to optimize the method.
+   *
+   * @param globalNodeNum  global node number of the node in the mesh _src_mesh
+   *
+   */
+  template<class MyMeshType>
+  inline void SplitterTetra<MyMeshType>::calculateNode2(typename MyMeshType::MyConnType globalNodeNum, const double* node)
+  {
+    double* transformedNode = new double[MyMeshType::MY_SPACEDIM];
+    assert(transformedNode != 0);
+    _t->apply(transformedNode, node);
+    _nodes[globalNodeNum] = transformedNode;
+  }
+
   /**
    * Calculates the volume contribution from the given TransformedTriangle and stores it with the given key in .
    * Calls TransformedTriangle::calculateIntersectionVolume to perform the calculation.
@@ -318,32 +521,50 @@ namespace INTERP_KERNEL
     _volumes.insert(std::make_pair(key, vol));
   }
 
+  /**
+   * Calculates the surface contribution from the given TransformedTriangle and stores it with the given key in.
+   * Calls TransformedTriangle::calculateIntersectionSurface to perform the calculation.
+   *
+   * @param tri    triangle for which to calculate the surface contribution
+   * @param key    key associated with the face
+   */
   template<class MyMeshType>
+  inline void SplitterTetra<MyMeshType>::calculateSurface(TransformedTriangle& tri, const TriangleFaceKey& key)
+  {
+    const double surf = tri.calculateIntersectionSurface(_t);
+    _volumes.insert(std::make_pair(key, surf));
+  }
+
+  template<class MyMeshTypeT, class MyMeshTypeS=MyMeshTypeT>
   class SplitterTetra2
   {
   public:
-    SplitterTetra2(const MyMeshType& targetMesh, const MyMeshType& srcMesh, SplittingPolicy policy);
+    SplitterTetra2(const MyMeshTypeT& targetMesh, const MyMeshTypeS& srcMesh, SplittingPolicy policy);
     ~SplitterTetra2();
     void releaseArrays();
-    void splitTargetCell(typename MyMeshType::MyConnType targetCell, typename MyMeshType::MyConnType nbOfNodesT,
-                         typename std::vector< SplitterTetra<MyMeshType>* >& tetra);
-    void fiveSplit(const int* const subZone, typename std::vector< SplitterTetra<MyMeshType>* >& tetra);
-    void sixSplit(const int* const subZone, typename std::vector< SplitterTetra<MyMeshType>* >& tetra);
-    void calculateGeneral24Tetra(typename std::vector< SplitterTetra<MyMeshType>* >& tetra);
-    void calculateGeneral48Tetra(typename std::vector< SplitterTetra<MyMeshType>* >& tetra);
-    void calculateSubNodes(const MyMeshType& targetMesh, typename MyMeshType::MyConnType targetCell);
-    inline const double* getCoordsOfSubNode(typename MyMeshType::MyConnType node);
-    inline const double* getCoordsOfSubNode2(typename MyMeshType::MyConnType node, typename MyMeshType::MyConnType& nodeId);
-    template<int n>
-    inline void calcBarycenter(double* barycenter, const typename MyMeshType::MyConnType* pts);
+    void splitTargetCell2(typename MyMeshTypeT::MyConnType targetCell, typename std::vector< SplitterTetra<MyMeshTypeS>* >& tetra);
+    void splitTargetCell(typename MyMeshTypeT::MyConnType targetCell, typename MyMeshTypeT::MyConnType nbOfNodesT,
+                         typename std::vector< SplitterTetra<MyMeshTypeS>* >& tetra);//to suppress
+    void fiveSplit(const int* const subZone, typename std::vector< SplitterTetra<MyMeshTypeS>* >& tetra);//to suppress
+    void sixSplit(const int* const subZone, typename std::vector< SplitterTetra<MyMeshTypeS>* >& tetra);//to suppress
+    void calculateGeneral24Tetra(typename std::vector< SplitterTetra<MyMeshTypeS>* >& tetra);//to suppress
+    void calculateGeneral48Tetra(typename std::vector< SplitterTetra<MyMeshTypeS>* >& tetra);//to suppress
+    void splitPyram5(typename std::vector< SplitterTetra<MyMeshTypeS>* >& tetra);//to suppress
+    void splitConvex(typename MyMeshTypeT::MyConnType                     targetCell,//to suppress
+                     typename std::vector< SplitterTetra<MyMeshTypeS>* >& tetra);//to suppress
+    void calculateSubNodes(const MyMeshTypeT& targetMesh, typename MyMeshTypeT::MyConnType targetCell);//to suppress
+    inline const double* getCoordsOfSubNode(typename MyMeshTypeT::MyConnType node);//to suppress
+    inline const double* getCoordsOfSubNode2(typename MyMeshTypeT::MyConnType node, typename MyMeshTypeT::MyConnType& nodeId);//to suppress
+    //template<int n>
+    inline void calcBarycenter(int n, double* barycenter, const typename MyMeshTypeT::MyConnType* pts);//to suppress
   private:
-    const MyMeshType& _target_mesh;
-    const MyMeshType& _src_mesh;
+    const MyMeshTypeT& _target_mesh;
+    const MyMeshTypeS& _src_mesh;
     SplittingPolicy _splitting_pol;
     /// vector of pointers to double[3] containing the coordinates of the
     /// (sub) - nodes of split target cell
     std::vector<const double*> _nodes;
-    std::vector<typename MyMeshType::MyConnType> _node_ids;
+    std::vector<typename MyMeshTypeT::MyConnType> _node_ids;
   };
 
   /**
@@ -353,9 +574,9 @@ namespace INTERP_KERNEL
    * @param  barycenter  pointer to double[3] array in which to store the result
    * @param  pts pointer to int[n] array containing the (sub)-nodes for which to calculate the barycenter
    */
-  template<class MyMeshType>
-  template<int n>
-  inline void SplitterTetra2<MyMeshType>::calcBarycenter(double* barycenter, const typename MyMeshType::MyConnType* pts)
+  template<class MyMeshTypeT, class MyMeshTypeS>
+  //template<int n>
+  inline void SplitterTetra2<MyMeshTypeT, MyMeshTypeS>::calcBarycenter(int n, double* barycenter, const typename MyMeshTypeT::MyConnType* pts)
   {
     barycenter[0] = barycenter[1] = barycenter[2] = 0.0;
     for(int i = 0; i < n ; ++i)
@@ -377,8 +598,8 @@ namespace INTERP_KERNEL
    * @param  node  local number of the (sub)-node 0,..,7 are the elements nodes, sub-nodes are numbered from 8,..
    * @return pointer to double[3] containing the coordinates of the nodes
    */
-  template<class MyMeshType>
-  inline const double* SplitterTetra2<MyMeshType>::getCoordsOfSubNode(typename MyMeshType::MyConnType node)
+  template<class MyMeshTypeT, class MyMeshTypeS>
+  inline const double* SplitterTetra2<MyMeshTypeT, MyMeshTypeS>::getCoordsOfSubNode(typename MyMeshTypeT::MyConnType node)
   {
     // replace "at()" with [] for unsafe but faster access
     return _nodes.at(node);
@@ -391,8 +612,8 @@ namespace INTERP_KERNEL
    * @param nodeId is an output that is node id in target whole mesh in C mode.
    * @return pointer to double[3] containing the coordinates of the nodes
    */
-  template<class MyMeshType>
-  const double* SplitterTetra2<MyMeshType>::getCoordsOfSubNode2(typename MyMeshType::MyConnType node, typename MyMeshType::MyConnType& nodeId)
+  template<class MyMeshTypeT, class MyMeshTypeS>
+  const double* SplitterTetra2<MyMeshTypeT, MyMeshTypeS>::getCoordsOfSubNode2(typename MyMeshTypeT::MyConnType node, typename MyMeshTypeT::MyConnType& nodeId)
   {
     const double *ret=_nodes.at(node);
     if(node<8)