PR: double nodes and flat elements for ASTER calculations in progress

idl/SMESH_MeshEditor.idl

@@ -1001,6 +1001,19 @@ module SMESH
                               out SMESH_Mesh     mesh,
                               out SMESH_Group    group) raises (SALOME::SALOME_Exception);
 
+    /*!
+     * \brief Double nodes on shared faces between groups of volumes and create flat elements on demand.
+     * The list of groups must describe a partition of the mesh volumes.
+     * The nodes of the internal faces at the boundaries of the groups are doubled.
+     * In option, the internal faces are replaced by flat elements.
+     * Triangles are transformed in prisms, and quadrangles in hexahedrons.
+     * \param theDomains - list of groups of volumes
+     * \param createJointElems - if TRUE, create the elements
+     * \return TRUE if operation has been completed successfully, FALSE otherwise
+     */
+    boolean DoubleNodesOnGroupBoundaries( in ListOfGroups theDomains,
+                                          in boolean createJointElems );
+
   };
 };
 

src/SMDS/SMDS_Downward.cxx

@@ -434,7 +434,7 @@ int SMDS_Down1D::computeFaces(int* pts, int* vtkIds, int nbcells, int* downFaces
           downTypes[cnt] = vtkType;
           cnt++;
         }
-      else
+      else if (SMDS_Downward::getCellDimension(vtkType) == 3)
         {
           int volId = _grid->CellIdToDownId(vtkId);
           SMDS_Downward * downvol = _grid->getDownArray(vtkType);
@@ -965,15 +965,15 @@ void SMDS_DownQuadTriangle::addDownCell(int cellId, int lowCellId, unsigned char
 {
   //ASSERT((cellId >=0)&& (cellId < _maxId));
   //ASSERT(aType == VTK_QUADRATIC_EDGE);
-  int *faces = &_cellIds[_nbDownCells * cellId];
+  int *edges = &_cellIds[_nbDownCells * cellId];
   for (int i = 0; i < _nbDownCells; i++)
     {
-      if (faces[i] < 0)
+      if (edges[i] < 0)
         {
-          faces[i] = lowCellId;
+          edges[i] = lowCellId;
           return;
         }
-      if (faces[i] == lowCellId)
+      if (edges[i] == lowCellId)
         return;
     }
   ASSERT(0);
@@ -1228,7 +1228,34 @@ SMDS_DownQuadTetra::~SMDS_DownQuadTetra()
 
 void SMDS_DownQuadTetra::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
 {
-  // TODO
+  set<int> setNodes;
+  setNodes.clear();
+  for (int i = 0; i < orderedNodes.size(); i++)
+    setNodes.insert(orderedNodes[i]);
+  //MESSAGE("cellId = " << cellId);
+
+  vtkIdType npts = 0;
+  vtkIdType *nodes; // will refer to the point id's of the volume
+  _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
+
+  set<int> tofind;
+  int ids[24] = { 0, 1, 2, 4, 5, 6,  0, 3, 1, 7, 8, 4,  2, 3, 0, 9, 7, 6,  1, 3, 2, 8, 9, 5 };
+//int ids[24] = { 2, 1, 0, 5, 4, 6,  1, 3, 0, 8, 7, 4,  0, 3, 2, 7, 9, 6,  2, 3, 1, 9, 8, 5 };
+  for (int k = 0; k < 4; k++)
+    {
+      tofind.clear();
+      for (int i = 0; i < 6; i++)
+        tofind.insert(nodes[ids[6 * k + i]]);
+      if (setNodes == tofind)
+        {
+          for (int i = 0; i < 6; i++)
+            orderedNodes[i] = nodes[ids[6 * k + i]];
+          return;
+        }
+    }
+  MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
+  MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2]);
+  MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
 }
 
 void SMDS_DownQuadTetra::addDownCell(int cellId, int lowCellId, unsigned char aType)
@@ -1282,8 +1309,8 @@ void SMDS_DownQuadTetra::computeFacesWithNodes(int cellId, ListElemByNodesType&
   facesWithNodes.elems[1].nodeIds[1] = nodes[1];
   facesWithNodes.elems[1].nodeIds[2] = nodes[3];
   facesWithNodes.elems[1].nodeIds[3] = nodes[4];
-  facesWithNodes.elems[1].nodeIds[4] = nodes[7];
+  facesWithNodes.elems[1].nodeIds[4] = nodes[8];
-  facesWithNodes.elems[1].nodeIds[5] = nodes[8];
+  facesWithNodes.elems[1].nodeIds[5] = nodes[7];
   facesWithNodes.elems[1].nbNodes = 6;
   facesWithNodes.elems[1].vtkType = VTK_QUADRATIC_TRIANGLE;
 
@@ -1291,8 +1318,8 @@ void SMDS_DownQuadTetra::computeFacesWithNodes(int cellId, ListElemByNodesType&
   facesWithNodes.elems[2].nodeIds[1] = nodes[2];
   facesWithNodes.elems[2].nodeIds[2] = nodes[3];
   facesWithNodes.elems[2].nodeIds[3] = nodes[6];
-  facesWithNodes.elems[2].nodeIds[4] = nodes[7];
+  facesWithNodes.elems[2].nodeIds[4] = nodes[9];
-  facesWithNodes.elems[2].nodeIds[5] = nodes[9];
+  facesWithNodes.elems[2].nodeIds[5] = nodes[7];
   facesWithNodes.elems[2].nbNodes = 6;
   facesWithNodes.elems[2].vtkType = VTK_QUADRATIC_TRIANGLE;
 
@@ -1300,8 +1327,8 @@ void SMDS_DownQuadTetra::computeFacesWithNodes(int cellId, ListElemByNodesType&
   facesWithNodes.elems[3].nodeIds[1] = nodes[2];
   facesWithNodes.elems[3].nodeIds[2] = nodes[3];
   facesWithNodes.elems[3].nodeIds[3] = nodes[5];
-  facesWithNodes.elems[3].nodeIds[4] = nodes[8];
+  facesWithNodes.elems[3].nodeIds[4] = nodes[9];
-  facesWithNodes.elems[3].nodeIds[5] = nodes[9];
+  facesWithNodes.elems[3].nodeIds[5] = nodes[8];
   facesWithNodes.elems[3].nbNodes = 6;
   facesWithNodes.elems[3].vtkType = VTK_QUADRATIC_TRIANGLE;
 }
@@ -1916,7 +1943,35 @@ SMDS_DownQuadHexa::~SMDS_DownQuadHexa()
 
 void SMDS_DownQuadHexa::getOrderedNodesOfFace(int cellId, std::vector<vtkIdType>& orderedNodes)
 {
-  // TODO
+  set<int> setNodes;
+  setNodes.clear();
+  for (int i = 0; i < orderedNodes.size(); i++)
+    setNodes.insert(orderedNodes[i]);
+  //MESSAGE("cellId = " << cellId);
+
+  vtkIdType npts = 0;
+  vtkIdType *nodes; // will refer to the point id's of the volume
+  _grid->GetCellPoints(this->_vtkCellIds[cellId], npts, nodes);
+
+  set<int> tofind;
+  //int ids[24] = { 3, 2, 1, 0,  4, 5, 6, 7,  7, 3, 0, 4,  4, 0, 1, 5,  5, 1, 2, 6,  6, 2, 3, 7};
+  int ids[48] = { 3, 2, 1, 0,10, 9, 8,11,   4, 5, 6, 7,12,13,14,15,   7, 3, 0, 4,19,11,16,15,
+                  4, 0, 1, 5,16, 8,17,12,   5, 1, 2, 6,17, 9,18,13,   6, 2, 3, 7,18,10,19,14};
+  for (int k = 0; k < 6; k++)
+    {
+      tofind.clear();
+      for (int i = 0; i < 8; i++)
+        tofind.insert(nodes[ids[8 * k + i]]);
+      if (setNodes == tofind)
+        {
+          for (int i = 0; i < 8; i++)
+            orderedNodes[i] = nodes[ids[8 * k + i]];
+          return;
+        }
+    }
+  MESSAGE("=== Problem volume " << _vtkCellIds[cellId] << " " << _grid->_mesh->fromVtkToSmds(_vtkCellIds[cellId]));
+  MESSAGE(orderedNodes[0] << " " << orderedNodes[1] << " " << orderedNodes[2] << " " << orderedNodes[3]);
+  MESSAGE(nodes[0] << " " << nodes[1] << " " << nodes[2] << " " << nodes[3]);
 }
 
 void SMDS_DownQuadHexa::addDownCell(int cellId, int lowCellId, unsigned char aType)

src/SMDS/SMDS_UnstructuredGrid.cxx

@@ -12,6 +12,7 @@
 #include <vtkUnsignedCharArray.h>
 
 #include <list>
+#include <climits>
 
 using namespace std;
 
@@ -805,3 +806,94 @@ void SMDS_UnstructuredGrid::BuildLinks()
   this->Links->BuildLinks(this, this->Connectivity);
   this->Links->Delete();
 }
+
+/*! Create a volume (prism or hexahedron) by duplication of a face.
+ * Designed for use in creation of flat elements separating volume domains.
+ * A face separating two domains is shared by two volume cells.
+ * All the nodes are already created (for the two faces).
+ * Each original Node is associated to corresponding nodes in the domains.
+ * Some nodes may be duplicated for more than two domains, when domain separations intersect.
+ * In that case, even some of the nodes to use for the original face may be changed.
+ * @param vtkVolId: vtk id of a volume containing the face, to get an orientation for the face.
+ * @param domain1: domain of the original face
+ * @param domain2: domain of the duplicated face
+ * @param originalNodes: the vtk node ids of the original face
+ * @param nodeDomains: map(original id --> map(domain --> duplicated node id))
+ * @return ok if success.
+ */
+bool SMDS_UnstructuredGrid::extrudeVolumeFromFace(int vtkVolId,
+                                                  int domain1,
+                                                  int domain2,
+                                                  std::set<int>& originalNodes,
+                                                  std::map<int, std::map<int, int> >& nodeDomains,
+                                                  std::map<int, std::map<long, int> >& nodeQuadDomains)
+{
+  //MESSAGE("extrudeVolumeFromFace " << vtkVolId);
+  vector<vtkIdType> orderedOriginals;
+  orderedOriginals.clear();
+  set<int>::const_iterator it = originalNodes.begin();
+  for (; it != originalNodes.end(); ++it)
+    orderedOriginals.push_back(*it);
+
+  int nbNodes = this->getOrderedNodesOfFace(vtkVolId, orderedOriginals);
+  vector<vtkIdType> orderedNodes;
+
+  switch (orderedOriginals.size())
+  {
+    case 3:
+    case 4:
+      for (int i = 0; i < nbNodes; i++)
+        orderedNodes.push_back(nodeDomains[orderedOriginals[i]][domain1]);
+      for (int i = 0; i < nbNodes; i++)
+        orderedNodes.push_back(nodeDomains[orderedOriginals[i]][domain2]);
+      break;
+    case 6:
+    case 8:
+      {
+        long dom1 = domain1;
+        long dom2 = domain2;
+        long dom1_2; // for nodeQuadDomains
+        if (domain1 < domain2)
+          dom1_2 = dom1 + INT_MAX * dom2;
+        else
+          dom1_2 = dom2 + INT_MAX * dom1;
+        //cerr << "dom1=" << dom1 << " dom2=" << dom2 << " dom1_2=" << dom1_2 << endl;
+        int ima = orderedOriginals.size();
+        int mid = orderedOriginals.size() / 2;
+        //cerr << "ima=" << ima << " mid=" << mid << endl;
+        for (int i = 0; i < mid; i++)
+          orderedNodes.push_back(nodeDomains[orderedOriginals[i]][domain1]);
+        for (int i = 0; i < mid; i++)
+          orderedNodes.push_back(nodeDomains[orderedOriginals[i]][domain2]);
+        for (int i = mid; i < ima; i++)
+          orderedNodes.push_back(nodeDomains[orderedOriginals[i]][domain1]);
+        for (int i = mid; i < ima; i++)
+          orderedNodes.push_back(nodeDomains[orderedOriginals[i]][domain2]);
+        for (int i = 0; i < mid; i++)
+          {
+            int oldId = orderedOriginals[i];
+            int newId;
+            if (nodeQuadDomains.count(oldId) && nodeQuadDomains[oldId].count(dom1_2))
+              newId = nodeQuadDomains[oldId][dom1_2];
+            else
+              {
+                double *coords = this->GetPoint(oldId);
+                SMDS_MeshNode *newNode = _mesh->AddNode(coords[0], coords[1], coords[2]);
+                newId = newNode->getVtkId();
+                std::map<long, int> emptyMap;
+                nodeQuadDomains[oldId] = emptyMap;
+                nodeQuadDomains[oldId][dom1_2] = newId;
+              }
+            orderedNodes.push_back(newId);
+          }
+      }
+      break;
+    default:
+      ASSERT(0);
+  }
+
+  SMDS_MeshVolume *vol = _mesh->AddVolumeFromVtkIds(orderedNodes);
+
+  // TODO update subshape list of elements and nodes
+  return vol;
+}

src/SMDS/SMDS_UnstructuredGrid.hxx

@@ -66,6 +66,12 @@ public:
   void ModifyCellNodes(int vtkVolId, std::map<int, int> localClonedNodeIds);
   int getOrderedNodesOfFace(int vtkVolId, std::vector<vtkIdType>& orderedNodes);
   void BuildLinks();
+  bool extrudeVolumeFromFace(int vtkVolId,
+                             int domain1,
+                             int domain2,
+                             std::set<int>& originalNodes,
+                             std::map<int,std::map<int,int> >& nodeDomains,
+                             std::map<int,std::map<long,int> >& nodeQuadDomains);
   vtkCellLinks* GetLinks()
   {
     return Links;

src/SMESH/SMESH_MeshEditor.cxx

@@ -10774,6 +10774,10 @@ bool SMESH_MeshEditor::DoubleNodesOnGroupBoundaries( const std::vector<TIDSorted
   //     get node id's of the face (id SMDS = id VTK)
   //     create flat element with old and new nodes if requested
 
+  // --- new quad nodes on flat quad elements: oldId --> ((domain1 X domain2) --> newId)
+  //     (domain1 X domain2) = domain1 + MAXINT*domain2
+  std::map<int, std::map<long,int> > nodeQuadDomains;
+
   if (createJointElems)
     {
       itface = faceDomains.begin();
@@ -10791,7 +10795,7 @@ bool SMESH_MeshEditor::DoubleNodesOnGroupBoundaries( const std::vector<TIDSorted
           int vtkVolId = itdom->second;
           itdom++;
           int dom2 = itdom->first;
-          meshDS->extrudeVolumeFromFace(vtkVolId, dom1, dom2, oldNodes, nodeDomains);
+          grid->extrudeVolumeFromFace(vtkVolId, dom1, dom2, oldNodes, nodeDomains, nodeQuadDomains);
         }
     }
 

src/SMESHDS/SMESHDS_Mesh.cxx

@@ -1975,42 +1975,3 @@ bool SMESHDS_Mesh::ModifyCellNodes(int vtkVolId, std::map<int,int> localClonedNo
   myGrid->ModifyCellNodes(vtkVolId, localClonedNodeIds);
   return true;
 }
-
-/*! Create a volume (prism or hexahedron) by duplication of a face.
- * Designed for use in creation of flat elements separating volume domains.
- * A face separating two domains is shared by two volume cells.
- * All the nodes are already created (for the two faces).
- * Each original Node is associated to corresponding nodes in the domains.
- * Some nodes may be duplicated for more than two domains, when domain separations intersect.
- * In that case, even some of the nodes to use for the original face may be changed.
- * @param vtkVolId: vtk id of a volume containing the face, to get an orientation for the face.
- * @param domain1: domain of the original face
- * @param domain2: domain of the duplicated face
- * @param originalNodes: the vtk node ids of the original face
- * @param nodeDomains: map(original id --> map(domain --> duplicated node id))
- * @return ok if success.
- */
-bool SMESHDS_Mesh::extrudeVolumeFromFace(int vtkVolId,
-                                         int domain1,
-                                         int domain2,
-                                         std::set<int>& originalNodes,
-                                         std::map<int,std::map<int,int> >& nodeDomains)
-{
-  //MESSAGE("extrudeVolumeFromFace " << vtkVolId);
-  vector<vtkIdType> orderedOriginals;
-  orderedOriginals.clear();
-  set<int>::const_iterator it = originalNodes.begin();
-  for (; it != originalNodes.end(); ++it)
-    orderedOriginals.push_back(*it);
-
-  int nbNodes = myGrid->getOrderedNodesOfFace(vtkVolId, orderedOriginals);
-  vector<vtkIdType> orderedNodes;
-  for (int i=0; i<nbNodes; i++)
-    orderedNodes.push_back(nodeDomains[orderedOriginals[i]][domain1]);
-  for (int i=0; i<nbNodes; i++)
-    orderedNodes.push_back(nodeDomains[orderedOriginals[i]][domain2]);
-  SMDS_MeshVolume *vol = this->AddVolumeFromVtkIds(orderedNodes);
-
-  // TODO update subshape list of elements and nodes
-  return vol;
-}

src/SMESHDS/SMESHDS_Mesh.hxx

@@ -401,11 +401,6 @@ public:
                              std::vector<const SMDS_MeshNode*> nodes,
                              std::vector<int>                  quantities);
   bool ModifyCellNodes(int smdsVolId, std::map<int,int> localClonedNodeIds);
-  bool extrudeVolumeFromFace(int vtkVolId,
-                             int domain1,
-                             int domain2,
-                             std::set<int>& originalNodes,
-                             std::map<int,std::map<int,int> >& nodeDomains);
   void Renumber (const bool isNodes, const int startID=1, const int deltaID=1);
 
   void SetNodeInVolume(SMDS_MeshNode * aNode, const TopoDS_Shell & S);