From 97d54079149b45691157ef0be3428e6b0144bafc Mon Sep 17 00:00:00 2001 From: ageay Date: Fri, 11 Jan 2013 10:49:07 +0000 Subject: [PATCH] MEDCouplingUMesh::findAndCorrectBadOriented3DCells --- src/MEDCoupling/MEDCouplingUMesh.cxx | 221 +++++++++++++++++------ src/MEDCoupling/MEDCouplingUMesh.hxx | 5 + src/MEDCoupling_Swig/MEDCouplingCommon.i | 2 + 3 files changed, 174 insertions(+), 54 deletions(-) diff --git a/src/MEDCoupling/MEDCouplingUMesh.cxx b/src/MEDCoupling/MEDCouplingUMesh.cxx index 9502459d1..465419372 100644 --- a/src/MEDCoupling/MEDCouplingUMesh.cxx +++ b/src/MEDCoupling/MEDCouplingUMesh.cxx @@ -235,6 +235,13 @@ void MEDCouplingUMesh::insertNextCell(INTERP_KERNEL::NormalizedCellType type, in throw INTERP_KERNEL::Exception("MEDCouplingUMesh::insertNextCell : nodal connectivity not set ! invoke allocateCells before calling insertNextCell !"); if((int)cm.getDimension()==_mesh_dim) { + if(!cm.isDynamic()) + if(size!=(int)cm.getNumberOfNodes()) + { + std::ostringstream oss; oss << "MEDCouplingUMesh::insertNextCell : Trying to push a " << cm.getRepr() << " cell with a size of " << size; + oss << " ! Expecting " << cm.getNumberOfNodes() << " !"; + throw INTERP_KERNEL::Exception(oss.str().c_str()); + } int idx=_nodal_connec_index->back(); int val=idx+size+1; _nodal_connec_index->pushBackSilent(val); @@ -3840,54 +3847,6 @@ DataArrayInt *MEDCouplingUMesh::convexEnvelop2D() throw(INTERP_KERNEL::Exception return isChanged.retn(); } -/*! - * This method is expected to be applied on a mesh with spaceDim==3 and meshDim==3. If not an exception will be thrown. - * This method analyzes only linear extruded 3D cells (NORM_HEXA8,NORM_PENTA6,NORM_HEXGP12...) - * If some extruded cells does not fulfill the MED norm for extruded cells (first face of 3D cell should be oriented to the exterior of the 3D cell). - * Some viewers are very careful of that (SMESH), but ParaVis ignore that. - */ -DataArrayInt *MEDCouplingUMesh::findAndCorrectBadOriented3DExtrudedCells() throw(INTERP_KERNEL::Exception) -{ - const char msg[]="check3DCellsWellOriented detection works only for 3D cells !"; - if(getMeshDimension()!=3) - throw INTERP_KERNEL::Exception(msg); - int spaceDim=getSpaceDimension(); - if(spaceDim!=3) - throw INTERP_KERNEL::Exception(msg); - // - int nbOfCells=getNumberOfCells(); - int *conn=_nodal_connec->getPointer(); - const int *connI=_nodal_connec_index->getConstPointer(); - const double *coo=getCoords()->getConstPointer(); - double vec0[3],vec1[3]; - MEDCouplingAutoRefCountObjectPtr cells(DataArrayInt::New()); - for(int i=0;i tmp=new int[connI[i+1]-connI[i]-1]; - int nbOfNodes=cm.fillSonCellNodalConnectivity(0,conn+connI[i]+1,tmp); - INTERP_KERNEL::areaVectorOfPolygon(tmp,nbOfNodes,coo,vec0); - const double *pt0=coo+3*conn[connI[i]+1]; - const double *pt1=coo+3*conn[connI[i]+nbOfNodes+1]; - vec1[0]=pt0[0]-pt1[0]; vec1[1]=pt0[1]-pt1[1]; vec1[2]=pt0[2]-pt1[2]; - double dot=vec0[0]*vec1[0]+vec0[1]*vec1[1]+vec0[2]*vec1[2]; - if(dot<0) - { - cells->pushBackSilent(i); - std::copy(conn+connI[i]+1,conn+connI[i+1],(int *)tmp); - for(int j=1;j& cell /*! * This method tries to orient correctly polhedrons cells. - * @throw when 'this' is not a mesh with meshdim==3 and spacedim==3. An exception is also thrown when the attempt of reparation fails. + * + * \throw when 'this' is not a mesh with meshdim==3 and spacedim==3. An exception is also thrown when the attempt of reparation fails. + * \sa MEDCouplingUMesh::findAndCorrectBadOriented3DCells */ void MEDCouplingUMesh::orientCorrectlyPolyhedrons() throw(INTERP_KERNEL::Exception) { @@ -4802,20 +4763,120 @@ void MEDCouplingUMesh::orientCorrectlyPolyhedrons() throw(INTERP_KERNEL::Excepti int *conn=_nodal_connec->getPointer(); const int *connI=_nodal_connec_index->getConstPointer(); const double *coordsPtr=_coords->getConstPointer(); - bool isModified=false; for(int i=0;igetPointer(); + const int *connI=_nodal_connec_index->getConstPointer(); + const double *coo=getCoords()->getConstPointer(); + MEDCouplingAutoRefCountObjectPtr cells(DataArrayInt::New()); cells->alloc(0,1); + for(int i=0;ipushBackSilent(i); + } + } + } + return cells.retn(); +} + +/*! + * This method is a faster method to correct orientation of all 3D cells in \a this. + * This method works only if \a this is a 3D mesh, that is to say a mesh with mesh dimension 3 and a space dimension 3. + * This method makes the hypothesis that \a this a coherent that is to say MEDCouplingUMesh::checkCoherency2 should throw no exception. + * + * \ret a newly allocated int array with one components containing cell ids renumbered to fit the convention of MED (MED file and MEDCoupling) + * \sa MEDCouplingUMesh::orientCorrectlyPolyhedrons, + */ +DataArrayInt *MEDCouplingUMesh::findAndCorrectBadOriented3DCells() throw(INTERP_KERNEL::Exception) +{ + if(getMeshDimension()!=3 || getSpaceDimension()!=3) + throw INTERP_KERNEL::Exception("Invalid mesh to apply findAndCorrectBadOriented3DCells on it : must be meshDim==3 and spaceDim==3 !"); + int nbOfCells=getNumberOfCells(); + int *conn=_nodal_connec->getPointer(); + const int *connI=_nodal_connec_index->getConstPointer(); + const double *coordsPtr=_coords->getConstPointer(); + MEDCouplingAutoRefCountObjectPtr ret=DataArrayInt::New(); ret->alloc(0,1); + for(int i=0;ipushBackSilent(i); + } + break; + } + case INTERP_KERNEL::NORM_PYRA5: + { + if(!IsPyra5WellOriented(conn+connI[i]+1,conn+connI[i+1],coordsPtr)) + { + std::swap(*(conn+connI[i]+2),*(conn+connI[i]+4)); + ret->pushBackSilent(i); + } + break; + } + case INTERP_KERNEL::NORM_PENTA6: + case INTERP_KERNEL::NORM_HEXA8: + case INTERP_KERNEL::NORM_HEXGP12: + { + if(!Is3DExtrudedStaticCellWellOriented(conn+connI[i]+1,conn+connI[i+1],coordsPtr)) + { + CorrectExtrudedStaticCell(conn+connI[i]+1,conn+connI[i+1]); + ret->pushBackSilent(i); + } + break; + } + case INTERP_KERNEL::NORM_POLYHED: { - TryToCorrectPolyhedronOrientation(conn+connI[i]+1,conn+connI[i+1],coordsPtr); - isModified=true; + if(!IsPolyhedronWellOriented(conn+connI[i]+1,conn+connI[i+1],coordsPtr)) + { + TryToCorrectPolyhedronOrientation(conn+connI[i]+1,conn+connI[i+1],coordsPtr); + ret->pushBackSilent(i); + } + break; } + default: + throw INTERP_KERNEL::Exception("MEDCouplingUMesh::orientCorrectly3DCells : Your mesh contains type of cell not supported yet ! send mail to anthony.geay@cea.fr to add it !"); + } } - if(isModified) - _nodal_connec->declareAsNew(); updateTime(); + return ret.retn(); } /*! @@ -6236,6 +6297,58 @@ bool MEDCouplingUMesh::IsPolyhedronWellOriented(const int *begin, const int *end return INTERP_KERNEL::calculateVolumeForPolyh2(begin,(int)std::distance(begin,end),coords)>-EPS_FOR_POLYH_ORIENTATION; } +/*! + * The 3D extruded static cell (PENTA6,HEXA8,HEXAGP12...) its connectivity nodes in [begin,end). + */ +bool MEDCouplingUMesh::Is3DExtrudedStaticCellWellOriented(const int *begin, const int *end, const double *coords) +{ + double vec0[3],vec1[3]; + std::size_t sz=std::distance(begin,end); + if(sz%2!=0) + throw INTERP_KERNEL::Exception("MEDCouplingUMesh::Is3DExtrudedStaticCellWellOriented : the length of nodal connectivity of extruded cell is not even !"); + int nbOfNodes=(int)sz/2; + INTERP_KERNEL::areaVectorOfPolygon(begin,nbOfNodes,coords,vec0); + const double *pt0=coords+3*begin[0]; + const double *pt1=coords+3*begin[nbOfNodes]; + vec1[0]=pt1[0]-pt0[0]; vec1[1]=pt1[1]-pt0[1]; vec1[2]=pt1[2]-pt0[2]; + return (vec0[0]*vec1[0]+vec0[1]*vec1[1]+vec0[2]*vec1[2])<0.; +} + +void MEDCouplingUMesh::CorrectExtrudedStaticCell(int *begin, int *end) +{ + std::size_t sz=std::distance(begin,end); + INTERP_KERNEL::AutoPtr tmp=new int[sz]; + std::size_t nbOfNodes(sz/2); + std::copy(begin,end,(int *)tmp); + for(std::size_t j=1;j(begin,4,coords,vec0); + const double *pt0=coords+3*begin[0],*pt1=coords+3*begin[4]; + return (vec0[0]*(pt1[0]-pt0[0])+vec0[1]*(pt1[1]-pt0[1])+vec0[2]*(pt1[2]-pt0[2]))<0.; +} + /*! * This method performs a simplyfication of a single polyedron cell. To do that each face of cell whose connectivity is defined by [\b begin, \b end) * is compared with the others in order to find faces in the same plane (with approx of eps). If any, the cells are grouped together and projected to diff --git a/src/MEDCoupling/MEDCouplingUMesh.hxx b/src/MEDCoupling/MEDCouplingUMesh.hxx index 4cfdfec52..534c7c7bc 100644 --- a/src/MEDCoupling/MEDCouplingUMesh.hxx +++ b/src/MEDCoupling/MEDCouplingUMesh.hxx @@ -164,6 +164,7 @@ namespace ParaMEDMEM MEDCOUPLING_EXPORT void checkButterflyCells(std::vector& cells, double eps=1e-12) const; MEDCOUPLING_EXPORT DataArrayInt *convexEnvelop2D() throw(INTERP_KERNEL::Exception); MEDCOUPLING_EXPORT DataArrayInt *findAndCorrectBadOriented3DExtrudedCells() throw(INTERP_KERNEL::Exception); + MEDCOUPLING_EXPORT DataArrayInt *findAndCorrectBadOriented3DCells() throw(INTERP_KERNEL::Exception); MEDCOUPLING_EXPORT void getBoundingBoxForBBTree(std::vector& bbox) const; MEDCOUPLING_EXPORT MEDCouplingUMesh *buildExtrudedMesh(const MEDCouplingUMesh *mesh1D, int policy); MEDCOUPLING_EXPORT bool isFullyQuadratic() const; @@ -217,6 +218,10 @@ namespace ParaMEDMEM MEDCOUPLING_EXPORT static void MergeNodesOnUMeshesSharingSameCoords(const std::vector& meshes, double eps) throw(INTERP_KERNEL::Exception); MEDCOUPLING_EXPORT static bool IsPolygonWellOriented(bool isQuadratic, const double *vec, const int *begin, const int *end, const double *coords); MEDCOUPLING_EXPORT static bool IsPolyhedronWellOriented(const int *begin, const int *end, const double *coords); + MEDCOUPLING_EXPORT static bool Is3DExtrudedStaticCellWellOriented(const int *begin, const int *end, const double *coords); + MEDCOUPLING_EXPORT static void CorrectExtrudedStaticCell(int *begin, int *end); + MEDCOUPLING_EXPORT static bool IsTetra4WellOriented(const int *begin, const int *end, const double *coords); + MEDCOUPLING_EXPORT static bool IsPyra5WellOriented(const int *begin, const int *end, const double *coords); MEDCOUPLING_EXPORT static void SimplifyPolyhedronCell(double eps, const DataArrayDouble *coords, const int *begin, const int *end, DataArrayInt *res) throw(INTERP_KERNEL::Exception); MEDCOUPLING_EXPORT static void ComputeVecAndPtOfFace(double eps, const double *coords, const int *begin, const int *end, double *v, double *p) throw(INTERP_KERNEL::Exception); MEDCOUPLING_EXPORT static void TryToCorrectPolyhedronOrientation(int *begin, int *end, const double *coords) throw(INTERP_KERNEL::Exception); diff --git a/src/MEDCoupling_Swig/MEDCouplingCommon.i b/src/MEDCoupling_Swig/MEDCouplingCommon.i index 9a30b57f6..41c1e107e 100644 --- a/src/MEDCoupling_Swig/MEDCouplingCommon.i +++ b/src/MEDCoupling_Swig/MEDCouplingCommon.i @@ -298,6 +298,7 @@ using namespace INTERP_KERNEL; %newobject ParaMEDMEM::MEDCouplingUMesh::keepCellIdsByType; %newobject ParaMEDMEM::MEDCouplingUMesh::Build0DMeshFromCoords; %newobject ParaMEDMEM::MEDCouplingUMesh::findAndCorrectBadOriented3DExtrudedCells; +%newobject ParaMEDMEM::MEDCouplingUMesh::findAndCorrectBadOriented3DCells; %newobject ParaMEDMEM::MEDCouplingUMesh::findCellIdsOnBoundary; %newobject ParaMEDMEM::MEDCouplingUMesh::computeSkin; %newobject ParaMEDMEM::MEDCouplingUMesh::getCellIdsLyingOnNodes; @@ -1412,6 +1413,7 @@ namespace ParaMEDMEM DataArrayInt *convexEnvelop2D() throw(INTERP_KERNEL::Exception); std::string cppRepr() const throw(INTERP_KERNEL::Exception); DataArrayInt *findAndCorrectBadOriented3DExtrudedCells() throw(INTERP_KERNEL::Exception); + DataArrayInt *findAndCorrectBadOriented3DCells() throw(INTERP_KERNEL::Exception); static MEDCouplingUMesh *Build0DMeshFromCoords(DataArrayDouble *da) throw(INTERP_KERNEL::Exception); static MEDCouplingUMesh *MergeUMeshes(const MEDCouplingUMesh *mesh1, const MEDCouplingUMesh *mesh2) throw(INTERP_KERNEL::Exception); static MEDCouplingUMesh *MergeUMeshesOnSameCoords(const MEDCouplingUMesh *mesh1, const MEDCouplingUMesh *mesh2) throw(INTERP_KERNEL::Exception); -- 2.39.2