]> SALOME platform Git repositories - tools/medcoupling.git/blobdiff - src/MEDCoupling/MEDCouplingUMesh.cxx
Salome HOME
Generalization of unstructured grid supported by the remapper.
[tools/medcoupling.git] / src / MEDCoupling / MEDCouplingUMesh.cxx
index fca3d1aa2a65063a3c2487722807e5616d3fde70..d594c82ab143ee3afe9f3fe8cdb16bf1dce804e0 100644 (file)
@@ -5239,6 +5239,7 @@ void MEDCouplingUMesh::tessellate2DCurve(double eps) throw(INTERP_KERNEL::Except
  *          and \a this->getMeshDimension() != 3. 
  *  \throw If \a policy is not one of the four discussed above.
  *  \throw If the nodal connectivity of cells is not defined.
+ * \sa MEDCouplingUMesh::tetrahedrize
  */
 DataArrayInt *MEDCouplingUMesh::simplexize(int policy) throw(INTERP_KERNEL::Exception)
 {
@@ -9240,7 +9241,6 @@ MEDCouplingUMesh *MEDCouplingUMesh::buildSpreadZonesWithPoly() const throw(INTER
  */
 std::vector<DataArrayInt *> MEDCouplingUMesh::partitionBySpreadZone() const throw(INTERP_KERNEL::Exception)
 {
-  //#if 0
   int nbOfCellsCur=getNumberOfCells();
   std::vector<DataArrayInt *> ret;
   if(nbOfCellsCur<=0)
@@ -9260,36 +9260,6 @@ std::vector<DataArrayInt *> MEDCouplingUMesh::partitionBySpreadZone() const thro
   for(std::vector< MEDCouplingAutoRefCountObjectPtr<DataArrayInt> >::iterator it=ret2.begin();it!=ret2.end();it++)
     ret.push_back((*it).retn());
   return ret;
-  //#endif
-#if 0
-  int nbOfCellsCur=getNumberOfCells();
-  DataArrayInt *neigh=0,*neighI=0;
-  computeNeighborsOfCells(neigh,neighI);
-  MEDCouplingAutoRefCountObjectPtr<DataArrayInt> neighAuto(neigh),neighIAuto(neighI);
-  MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ids=DataArrayInt::New(); ids->alloc(nbOfCellsCur,1); ids->iota();
-  std::vector<DataArrayInt *> ret;
-  std::vector< MEDCouplingAutoRefCountObjectPtr<DataArrayInt> > ret2;
-  while(nbOfCellsCur>0)
-    {
-      MEDCouplingAutoRefCountObjectPtr<DataArrayInt> tmp=MEDCouplingUMesh::ComputeSpreadZoneGradually(neighAuto,neighIAuto);
-      MEDCouplingAutoRefCountObjectPtr<DataArrayInt> tmp3=tmp->buildComplement(nbOfCellsCur);
-      MEDCouplingAutoRefCountObjectPtr<DataArrayInt> tmp2=ids->selectByTupleId(tmp->begin(),tmp->end());
-      ret2.push_back(tmp2);  ret.push_back(tmp2);
-      nbOfCellsCur=tmp3->getNumberOfTuples();
-      if(nbOfCellsCur>0)
-        {
-          ids=ids->selectByTupleId(tmp3->begin(),tmp3->end());
-          MEDCouplingUMesh::ExtractFromIndexedArrays(tmp3->begin(),tmp3->end(),neighAuto,neighIAuto,neigh,neighI);
-          neighAuto=neigh;
-          neighIAuto=neighI;
-          MEDCouplingAutoRefCountObjectPtr<DataArrayInt> renum=tmp3->invertArrayN2O2O2N(nbOfCellsCur+tmp->getNumberOfTuples());
-          neighAuto->transformWithIndArr(renum->begin(),renum->end());
-        }
-    }
-  for(std::vector<DataArrayInt *>::const_iterator it=ret.begin();it!=ret.end();it++)
-    (*it)->incrRef();
-  return ret;
-#endif
 }
 
 /*!
@@ -9316,6 +9286,78 @@ DataArrayInt *MEDCouplingUMesh::ComputeRangesFromTypeDistribution(const std::vec
   return ret.retn();
 }
 
+/*!
+ * This method expects that \a this a 3D mesh (spaceDim=3 and meshDim=3) with all coordinates and connectivities set.
+ * All cells in \a this are expected to be linear 3D cells.
+ * This method will split **all** 3D cells in \a this into INTERP_KERNEL::NORM_TETRA4 cells and put them in the returned mesh.
+ * It leads to an increase to number of cells.
+ * This method contrary to MEDCouplingUMesh::simplexize can append coordinates in \a this to perform its work.
+ * The \a nbOfAdditionalPoints returned value informs about it. If > 0, the coordinates array in returned mesh will have \a nbOfAdditionalPoints 
+ * more tuples (nodes) than in \a this. Anyway, all the nodes in \a this (with the same order) will be in the returned mesh.
+ *
+ * \param [in] policy - the policy of splitting that must be in (PLANAR_FACE_5, PLANAR_FACE_6, GENERAL_24, GENERAL_48). The policy will be used only for INTERP_KERNEL::NORM_HEXA8 cells.
+ *                      For all other cells, the splitting policy will be ignored.
+ * \param [out] nbOfAdditionalPoints - number of nodes added to \c this->_coords. If > 0 a new coordinates object will be constructed result of the aggregation of the old one and the new points added. 
+ * \param [out] n2oCells - A new instance of DataArrayInt holding, for each new cell,
+ *          an id of old cell producing it. The caller is to delete this array using
+ *         decrRef() as it is no more needed.
+ * \return MEDCoupling1SGTUMesh * - the mesh containing only INTERP_KERNEL::NORM_TETRA4 cells.
+ * 
+ * \throw If \a this is not a 3D mesh (spaceDim==3 and meshDim==3).
+ * \throw If \a this is not fully constituted with linear 3D cells.
+ * \sa MEDCouplingUMesh::simplexize
+ */
+MEDCoupling1SGTUMesh *MEDCouplingUMesh::tetrahedrize(int policy, DataArrayInt *& n2oCells, int& nbOfAdditionalPoints) const throw(INTERP_KERNEL::Exception)
+{
+  INTERP_KERNEL::SplittingPolicy pol((INTERP_KERNEL::SplittingPolicy)policy);
+  checkConnectivityFullyDefined();
+  if(getMeshDimension()!=3 || getSpaceDimension()!=3)
+    throw INTERP_KERNEL::Exception("MEDCouplingUMesh::tetrahedrize : only available for mesh with meshdim == 3 and spacedim == 3 !");
+  int nbOfCells(getNumberOfCells()),nbNodes(getNumberOfNodes());
+  MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> ret0(MEDCoupling1SGTUMesh::New(getName().c_str(),INTERP_KERNEL::NORM_TETRA4));
+  MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret(DataArrayInt::New()); ret->alloc(nbOfCells,1);
+  int *retPt(ret->getPointer());
+  MEDCouplingAutoRefCountObjectPtr<DataArrayInt> newConn(DataArrayInt::New()); newConn->alloc(0,1);
+  MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> addPts(DataArrayDouble::New()); addPts->alloc(0,1);
+  const int *oldc(_nodal_connec->begin());
+  const int *oldci(_nodal_connec_index->begin());
+  const double *coords(_coords->begin());
+  for(int i=0;i<nbOfCells;i++,oldci++,retPt++)
+    {
+      std::vector<int> a; std::vector<double> b;
+      INTERP_KERNEL::SplitIntoTetras(pol,(INTERP_KERNEL::NormalizedCellType)oldc[oldci[0]],oldc+oldci[0]+1,oldc+oldci[1],coords,a,b);
+      std::size_t nbOfTet(a.size()/4); *retPt=(int)nbOfTet;
+      const int *aa(&a[0]);
+      if(!b.empty())
+        {
+          for(std::vector<int>::iterator it=a.begin();it!=a.end();it++)
+            if(*it<0)
+              *it=(-(*(it))-1+nbNodes);
+          addPts->insertAtTheEnd(b.begin(),b.end());
+          nbNodes+=(int)b.size()/3;
+        }
+      for(std::size_t j=0;j<nbOfTet;j++,aa+=4)
+        newConn->insertAtTheEnd(aa,aa+4);
+    }
+  if(!addPts->empty())
+    {
+      addPts->rearrange(3);
+      nbOfAdditionalPoints=addPts->getNumberOfTuples();
+      addPts=DataArrayDouble::Aggregate(getCoords(),addPts);
+      ret0->setCoords(addPts);
+    }
+  else
+    {
+      nbOfAdditionalPoints=0;
+      ret0->setCoords(getCoords());
+    }
+  ret0->setNodalConnectivity(newConn);
+  //
+  ret->computeOffsets2();
+  n2oCells=ret->buildExplicitArrOfSliceOnScaledArr(0,nbOfCells,1);
+  return ret0.retn();
+}
+
 MEDCouplingUMeshCellIterator::MEDCouplingUMeshCellIterator(MEDCouplingUMesh *mesh):_mesh(mesh),_cell(new MEDCouplingUMeshCell(mesh)),
                                                                                    _own_cell(true),_cell_id(-1),_nb_cell(0)
 {