From: Anthony Geay <anthony.geay@edf.fr>
Date: Mon, 7 Aug 2023 15:46:11 +0000 (+0200)
Subject: [EDF27859] : GENERAL_24 management
X-Git-Url: http://git.salome-platform.org/gitweb/?a=commitdiff_plain;h=88e28cdfe7064e31b91238894fcc8f75cb861470;p=tools%2Fmedcoupling.git

[EDF27859] : GENERAL_24 management
---

diff --git a/src/INTERP_KERNEL/PolyhedronIntersectorP1P0.txx b/src/INTERP_KERNEL/PolyhedronIntersectorP1P0.txx
index 5741db4ac..7c4776d36 100644
--- a/src/INTERP_KERNEL/PolyhedronIntersectorP1P0.txx
+++ b/src/INTERP_KERNEL/PolyhedronIntersectorP1P0.txx
@@ -66,6 +66,23 @@ namespace INTERP_KERNEL
     _tetra.clear();
   }
 
+  template<class RowType, class ConnType>
+  void AddContributionInRow(RowType& row, ConnType colId, double value)
+  {
+    if(value != 0.)
+    {
+      typename RowType::const_iterator iterRes=row.find(colId);
+      if(iterRes==row.end())
+        row.insert(std::make_pair(colId,value));
+      else
+      {
+        double val=(*iterRes).second+value;
+        row.erase(colId);
+        row.insert(std::make_pair(colId,val));
+      }
+    }
+  }
+
   /**
    * Calculates the volume of intersection of an element in the source mesh and the target element
    * represented by the object.
@@ -80,37 +97,48 @@ namespace INTERP_KERNEL
   template<class MyMeshType, class MyMatrix>
   void PolyhedronIntersectorP1P0<MyMeshType,MyMatrix>::intersectCells(ConnType targetCell, const std::vector<ConnType>& srcCells, MyMatrix& res)
   {
-    SplitterTetra<MyMeshType>* subTetras[24];
     typename MyMatrix::value_type& resRow=res[targetCell];
-    for(typename std::vector<ConnType>::const_iterator iterCellS=srcCells.begin();iterCellS!=srcCells.end();iterCellS++)
-      {
-        releaseArrays();
-        ConnType nbOfNodesS=Intersector3D<MyMeshType,MyMatrix>::_src_mesh.getNumberOfNodesOfElement(OTT<ConnType,numPol>::indFC(*iterCellS));
-        _split.splitTargetCell(*iterCellS,nbOfNodesS,_tetra);
-        for(typename std::vector<SplitterTetra<MyMeshType>*>::iterator iter = _tetra.begin(); iter != _tetra.end(); ++iter)
-          {
-            (*iter)->splitIntoDualCells(subTetras);
-            for(int i=0;i<24;i++)
-              {
-                SplitterTetra<MyMeshType> *tmp=subTetras[i];
-                double volume = tmp->intersectSourceCell(targetCell);
-                ConnType sourceNode=tmp->getId(0);
-                if(volume!=0.)
-                  {
-                    typename MyMatrix::value_type::const_iterator iterRes=resRow.find(OTT<ConnType,numPol>::indFC(sourceNode));
-                    if(iterRes==resRow.end())
-                      resRow.insert(std::make_pair(OTT<ConnType,numPol>::indFC(sourceNode),volume));
-                    else
-                      {
-                        double val=(*iterRes).second+volume;
-                        resRow.erase(OTT<ConnType,numPol>::indFC(sourceNode));
-                        resRow.insert(std::make_pair(OTT<ConnType,numPol>::indFC(sourceNode),val));
-                      }
-                  }
-                delete tmp;
-              }
-          }
-      }
+    if( _split.getSplittingPolicy() != GENERAL_24)
+    {
+      SplitterTetra<MyMeshType>* subTetras[24];
+      for(typename std::vector<ConnType>::const_iterator iterCellS=srcCells.begin();iterCellS!=srcCells.end();iterCellS++)
+        {
+          releaseArrays();
+          ConnType nbOfNodesS=Intersector3D<MyMeshType,MyMatrix>::_src_mesh.getNumberOfNodesOfElement(OTT<ConnType,numPol>::indFC(*iterCellS));
+          _split.splitTargetCell(*iterCellS,nbOfNodesS,_tetra);
+          for(typename std::vector<SplitterTetra<MyMeshType>*>::const_iterator iter = _tetra.cbegin(); iter != _tetra.cend(); ++iter)
+            {
+              (*iter)->splitIntoDualCells(subTetras);
+              double vol2 = 0.;
+              for(int i=0;i<24;i++)
+                {
+                  SplitterTetra<MyMeshType> *tmp=subTetras[i];
+                  double volume = tmp->intersectSourceCell(targetCell);
+                  vol2 += volume;
+                  ConnType sourceNode=tmp->getId(0);
+                  AddContributionInRow(resRow,OTT<ConnType,numPol>::indFC(sourceNode),volume);
+                  delete tmp;
+                }
+            }
+        }
+    }
+    else
+    {
+      for(typename std::vector<ConnType>::const_iterator iterCellS=srcCells.begin();iterCellS!=srcCells.end();iterCellS++)
+        {
+          releaseArrays();
+          ConnType nbOfNodesS=Intersector3D<MyMeshType,MyMatrix>::_src_mesh.getNumberOfNodesOfElement(OTT<ConnType,numPol>::indFC(*iterCellS));
+          _split.splitTargetCell(*iterCellS,nbOfNodesS,_tetra);
+          for(typename std::vector<SplitterTetra<MyMeshType>*>::const_iterator iter = _tetra.cbegin(); iter != _tetra.cend(); ++iter)
+            {
+              double volume = std::abs( (*iter)->intersectSourceCell(targetCell) );
+              // node #0 is for internal node node #1 is for the node at the middle of the face
+              ConnType sourceNode0( (*iter)->getId(2) ), sourceNode1( (*iter)->getId(3) );
+              AddContributionInRow(resRow,OTT<ConnType,numPol>::indFC(sourceNode0),volume/2.);
+              AddContributionInRow(resRow,OTT<ConnType,numPol>::indFC(sourceNode1),volume/2.);
+            }
+        }
+    }
   }
 }
 
diff --git a/src/INTERP_KERNEL/SplitterTetra.hxx b/src/INTERP_KERNEL/SplitterTetra.hxx
index 30c46ddcb..b65bb4994 100644
--- a/src/INTERP_KERNEL/SplitterTetra.hxx
+++ b/src/INTERP_KERNEL/SplitterTetra.hxx
@@ -544,6 +544,7 @@ namespace INTERP_KERNEL
     SplitterTetra2(const MyMeshTypeT& targetMesh, const MyMeshTypeS& srcMesh, SplittingPolicy policy);
     ~SplitterTetra2();
     void releaseArrays();
+    SplittingPolicy getSplittingPolicy() const { return _splitting_pol; }
     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
diff --git a/src/INTERP_KERNEL/SplitterTetra.txx b/src/INTERP_KERNEL/SplitterTetra.txx
index 11202f04c..a0da9c3c0 100644
--- a/src/INTERP_KERNEL/SplitterTetra.txx
+++ b/src/INTERP_KERNEL/SplitterTetra.txx
@@ -1119,7 +1119,8 @@ namespace INTERP_KERNEL
     typename MyMeshTypeS::MyConnType conn[4];
     // get the cell center
     conn[0] = 14;
-    nodes[0] = getCoordsOfSubNode(conn[0]);
+    typename MyMeshTypeS::MyConnType realConn;
+    nodes[0] = getCoordsOfSubNode(conn[0]); 
 
     for(int faceCenterNode = 8; faceCenterNode < 14; ++faceCenterNode)
       {
@@ -1131,8 +1132,8 @@ namespace INTERP_KERNEL
             const int row = 4*(faceCenterNode - 8) + j;
             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]);
+            nodes[2] = getCoordsOfSubNode2(conn[2],realConn); conn[2] = realConn;
+            nodes[3] = getCoordsOfSubNode2(conn[3],realConn); conn[3] = realConn;
 
             SplitterTetra<MyMeshTypeS>* t = new SplitterTetra<MyMeshTypeS>(_src_mesh, nodes, conn);
             tetra.push_back(t);