From: ageay Date: Mon, 29 Nov 2010 12:04:56 +0000 (+0000) Subject: *** empty log message *** X-Git-Url: http://git.salome-platform.org/gitweb/?a=commitdiff_plain;h=0a9d87e3e5689552023b02fe391228549f7acfff;p=tools%2Fmedcoupling.git *** empty log message *** --- diff --git a/src/MEDCoupling/MEDCouplingMemArray.cxx b/src/MEDCoupling/MEDCouplingMemArray.cxx index 463012e9d..10e5f2256 100644 --- a/src/MEDCoupling/MEDCouplingMemArray.cxx +++ b/src/MEDCoupling/MEDCouplingMemArray.cxx @@ -156,6 +156,21 @@ void DataArrayDouble::fillWithValue(double val) declareAsNew(); } +bool DataArrayDouble::isUniform(double val, double eps) const +{ + if(getNumberOfComponents()!=1) + throw INTERP_KERNEL::Exception("DataArrayDouble::isUniform : must be applied on DataArrayDouble with only one component !"); + int nbOfTuples=getNumberOfTuples(); + const double *w=getConstPointer(); + const double *end=w+nbOfTuples; + const double vmin=val-eps; + const double vmax=val+eps; + for(;w!=end;w++) + if(*wvmax) + return false; + return true; +} + std::string DataArrayDouble::repr() const { std::ostringstream ret; @@ -227,6 +242,26 @@ DataArrayInt *DataArrayDouble::convertToIntArr() const return ret; } +DataArrayDouble *DataArrayDouble::fromNoInterlace() const throw(INTERP_KERNEL::Exception) +{ + if(_mem.isNull()) + throw INTERP_KERNEL::Exception("DataArrayDouble::fromNoInterlace : Not defined array !"); + double *tab=_mem.fromNoInterlace(getNumberOfComponents()); + DataArrayDouble *ret=DataArrayDouble::New(); + ret->useArray(tab,true,CPP_DEALLOC,getNumberOfTuples(),getNumberOfComponents()); + return ret; +} + +DataArrayDouble *DataArrayDouble::toNoInterlace() const throw(INTERP_KERNEL::Exception) +{ + if(_mem.isNull()) + throw INTERP_KERNEL::Exception("DataArrayDouble::fromNoInterlace : Not defined array !"); + double *tab=_mem.toNoInterlace(getNumberOfComponents()); + DataArrayDouble *ret=DataArrayDouble::New(); + ret->useArray(tab,true,CPP_DEALLOC,getNumberOfTuples(),getNumberOfComponents()); + return ret; +} + /*! * This method does \b not change the number of tuples after this call. * Only a permutation is done. If a permutation reduction is needed substr, or selectByTupleId should be used. @@ -537,6 +572,63 @@ double DataArrayDouble::accumulate(int compId) const return ret; } +DataArrayDouble *DataArrayDouble::fromPolarToCart() const +{ + int nbOfComp=getNumberOfComponents(); + if(nbOfComp!=2) + throw INTERP_KERNEL::Exception("DataArrayDouble::fromPolarToCart : must be an array with exactly 2 components !"); + int nbOfTuple=getNumberOfTuples(); + DataArrayDouble *ret=DataArrayDouble::New(); + ret->alloc(nbOfTuple,2); + double *w=ret->getPointer(); + const double *wIn=getConstPointer(); + for(int i=0;ialloc(getNumberOfTuples(),3); + double *w=ret->getPointer(); + const double *wIn=getConstPointer(); + for(int i=0;isetInfoOnComponent(2,getInfoOnComponent(2).c_str()); + return ret; +} + +DataArrayDouble *DataArrayDouble::fromSpherToCart() const +{ + int nbOfComp=getNumberOfComponents(); + if(nbOfComp!=3) + throw INTERP_KERNEL::Exception("DataArrayDouble::fromSpherToCart : must be an array with exactly 3 components !"); + int nbOfTuple=getNumberOfTuples(); + DataArrayDouble *ret=DataArrayDouble::New(); + ret->alloc(getNumberOfTuples(),3); + double *w=ret->getPointer(); + const double *wIn=getConstPointer(); + for(int i=0;iuseArray(tab,true,CPP_DEALLOC,getNumberOfTuples(),getNumberOfComponents()); + return ret; +} + +DataArrayInt *DataArrayInt::toNoInterlace() const throw(INTERP_KERNEL::Exception) +{ + if(_mem.isNull()) + throw INTERP_KERNEL::Exception("DataArrayInt::toNoInterlace : Not defined array !"); + int *tab=_mem.toNoInterlace(getNumberOfComponents()); + DataArrayInt *ret=DataArrayInt::New(); + ret->useArray(tab,true,CPP_DEALLOC,getNumberOfTuples(),getNumberOfComponents()); + return ret; +} + void DataArrayInt::renumberInPlace(const int *old2New) { int nbTuples=getNumberOfTuples(); @@ -1418,6 +1530,19 @@ bool DataArrayInt::isIdentity() const return true; } +bool DataArrayInt::isUniform(int val) const +{ + if(getNumberOfComponents()!=1) + throw INTERP_KERNEL::Exception("DataArrayInt::isUniform : must be applied on DataArrayInt with only one component !"); + int nbOfTuples=getNumberOfTuples(); + const int *w=getConstPointer(); + const int *end=w+nbOfTuples; + for(;w!=end;w++) + if(*w!=val) + return false; + return true; +} + DataArrayDouble *DataArrayInt::convertToDblArr() const { DataArrayDouble *ret=DataArrayDouble::New(); diff --git a/src/MEDCoupling/MEDCouplingMemArray.hxx b/src/MEDCoupling/MEDCouplingMemArray.hxx index c6cd1410b..2d0be0e17 100644 --- a/src/MEDCoupling/MEDCouplingMemArray.hxx +++ b/src/MEDCoupling/MEDCouplingMemArray.hxx @@ -54,6 +54,7 @@ namespace ParaMEDMEM public: MemArray():_nb_of_elem(-1),_ownership(false),_dealloc(CPP_DEALLOC) { } MemArray(const MemArray& other); + bool isNull() const { return _pointer.isNull(); } const T *getConstPointerLoc(int offset) const { return _pointer.getConstPointerLoc(offset); } const T *getConstPointer() const { return _pointer.getConstPointer(); } T *getPointer() const { return _pointer.getPointer(); } @@ -64,6 +65,8 @@ namespace ParaMEDMEM void repr(int sl, std::ostream& stream) const; void reprZip(int sl, std::ostream& stream) const; void fillWithValue(const T& val); + T *fromNoInterlace(int nbOfComp) const; + T *toNoInterlace(int nbOfComp) const; void alloc(int nbOfElements); void reAlloc(int newNbOfElements); void useArray(const T *array, bool ownership, DeallocType type, int nbOfElem); @@ -118,6 +121,7 @@ namespace ParaMEDMEM MEDCOUPLING_EXPORT void alloc(int nbOfTuple, int nbOfCompo); MEDCOUPLING_EXPORT void fillWithZero(); MEDCOUPLING_EXPORT void fillWithValue(double val); + MEDCOUPLING_EXPORT bool isUniform(double val, double eps) const; MEDCOUPLING_EXPORT std::string repr() const; MEDCOUPLING_EXPORT std::string reprZip() const; MEDCOUPLING_EXPORT void reprStream(std::ostream& stream) const; @@ -129,6 +133,8 @@ namespace ParaMEDMEM //!alloc or useArray should have been called before. MEDCOUPLING_EXPORT void reAlloc(int nbOfTuples); MEDCOUPLING_EXPORT DataArrayInt *convertToIntArr() const; + MEDCOUPLING_EXPORT DataArrayDouble *fromNoInterlace() const throw(INTERP_KERNEL::Exception); + MEDCOUPLING_EXPORT DataArrayDouble *toNoInterlace() const throw(INTERP_KERNEL::Exception); MEDCOUPLING_EXPORT void renumberInPlace(const int *old2New); MEDCOUPLING_EXPORT void renumberInPlaceR(const int *new2Old); MEDCOUPLING_EXPORT DataArrayDouble *renumber(const int *old2New) const; @@ -155,6 +161,9 @@ namespace ParaMEDMEM MEDCOUPLING_EXPORT double getAverageValue() const throw(INTERP_KERNEL::Exception); MEDCOUPLING_EXPORT void accumulate(double *res) const; MEDCOUPLING_EXPORT double accumulate(int compId) const; + MEDCOUPLING_EXPORT DataArrayDouble *fromPolarToCart() const; + MEDCOUPLING_EXPORT DataArrayDouble *fromCylToCart() const; + MEDCOUPLING_EXPORT DataArrayDouble *fromSpherToCart() const; MEDCOUPLING_EXPORT DataArrayDouble *doublyContractedProduct() const throw(INTERP_KERNEL::Exception); MEDCOUPLING_EXPORT DataArrayDouble *determinant() const throw(INTERP_KERNEL::Exception); MEDCOUPLING_EXPORT DataArrayDouble *eigenValues() const throw(INTERP_KERNEL::Exception); @@ -216,6 +225,8 @@ namespace ParaMEDMEM //!alloc or useArray should have been called before. MEDCOUPLING_EXPORT void reAlloc(int nbOfTuples); MEDCOUPLING_EXPORT DataArrayDouble *convertToDblArr() const; + MEDCOUPLING_EXPORT DataArrayInt *fromNoInterlace() const throw(INTERP_KERNEL::Exception); + MEDCOUPLING_EXPORT DataArrayInt *toNoInterlace() const throw(INTERP_KERNEL::Exception); MEDCOUPLING_EXPORT void renumberInPlace(const int *old2New); MEDCOUPLING_EXPORT void renumberInPlaceR(const int *new2Old); MEDCOUPLING_EXPORT DataArrayInt *renumber(const int *old2New) const; @@ -223,6 +234,7 @@ namespace ParaMEDMEM MEDCOUPLING_EXPORT DataArrayInt *renumberAndReduce(const int *old2NewBg, int newNbOfTuple) const; MEDCOUPLING_EXPORT DataArrayInt *selectByTupleId(const int *new2OldBg, const int *new2OldEnd) const; MEDCOUPLING_EXPORT bool isIdentity() const; + MEDCOUPLING_EXPORT bool isUniform(int val) const; MEDCOUPLING_EXPORT DataArrayInt *substr(int tupleIdBg, int tupleIdEnd=-1) const throw(INTERP_KERNEL::Exception); MEDCOUPLING_EXPORT DataArrayInt *changeNbOfComponents(int newNbOfComp, int dftValue) const throw(INTERP_KERNEL::Exception); MEDCOUPLING_EXPORT DataArrayInt *keepSelectedComponents(const std::vector& compoIds) const throw(INTERP_KERNEL::Exception); diff --git a/src/MEDCoupling/MEDCouplingMemArray.txx b/src/MEDCoupling/MEDCouplingMemArray.txx index 5a1a048d5..1c6e8cefc 100644 --- a/src/MEDCoupling/MEDCouplingMemArray.txx +++ b/src/MEDCoupling/MEDCouplingMemArray.txx @@ -180,6 +180,32 @@ namespace ParaMEDMEM T *pt=_pointer.getPointer(); std::fill(pt,pt+_nb_of_elem,val); } + + template + T *MemArray::fromNoInterlace(int nbOfComp) const + { + const T *pt=_pointer.getConstPointer(); + int nbOfTuples=_nb_of_elem/nbOfComp; + T *ret=new T[_nb_of_elem]; + T *w=ret; + for(int i=0;i + T *MemArray::toNoInterlace(int nbOfComp) const + { + const T *pt=_pointer.getConstPointer(); + int nbOfTuples=_nb_of_elem/nbOfComp; + T *ret=new T[_nb_of_elem]; + T *w=ret; + for(int i=0;i void MemArray::alloc(int nbOfElements) diff --git a/src/MEDCoupling/MEDCouplingUMesh.cxx b/src/MEDCoupling/MEDCouplingUMesh.cxx index 11f9d12ea..8b62afb68 100644 --- a/src/MEDCoupling/MEDCouplingUMesh.cxx +++ b/src/MEDCoupling/MEDCouplingUMesh.cxx @@ -563,6 +563,59 @@ void MEDCouplingUMesh::convertToPolyTypes(const std::vector& cellIdsToConve computeTypes(); } +/*! + * This method is the opposite of ParaMEDMEM::MEDCouplingUMesh::convertToPolyTypes method. + * The aim is to take all polygons or polyhedrons cell and to try to traduce them into classical cells. + * + */ +void MEDCouplingUMesh::unPolyze() +{ + checkFullyDefined(); + if(getMeshDimension()<=1) + throw INTERP_KERNEL::Exception("MEDCouplingUMesh::convertDegeneratedCells works on umeshes with meshdim equals to 2 or 3 !"); + int nbOfCells=getNumberOfCells(); + if(nbOfCells<1) + return ; + int initMeshLgth=getMeshLength(); + int *conn=_nodal_connec->getPointer(); + int *index=_nodal_connec_index->getPointer(); + int posOfCurCell=0; + int newPos=0; + int lgthOfCurCell; + for(int i=0;ireAlloc(newPos); + computeTypes(); +} + /*! * Array returned is the correspondance new to old. * The maximum value stored in returned array is the number of nodes of 'this' minus 1 after call of this method. @@ -1886,9 +1939,6 @@ void MEDCouplingUMesh::getCellsContainingPoints(const double *pos, int nbOfPoint else throw INTERP_KERNEL::Exception("For spaceDim==2 only meshDim==2 implemented for getelementscontainingpoints !"); } - - - } /*! @@ -2156,6 +2206,46 @@ void MEDCouplingUMesh::convertQuadraticCellsToLinear() throw(INTERP_KERNEL::Exce setConnectivity(newConn,newConnI,false); } +/*! + * This method converts all degenerated cells to simpler cells. For example a NORM_QUAD4 cell consituted from 2 same node id in its + * nodal connectivity will be transform to a NORM_TRI3 cell. + * This method works \b only \b on \b linear cells. + * This method works on nodes ids, that is to say a call to ParaMEDMEM::MEDCouplingUMesh::mergeNodes + * method could be usefull before calling this method in case of presence of several pair of nodes located on same position. + * This method throws an exception if 'this' is not fully defined (connectivity). + * This method throws an exception too if a "too" degenerated cell is detected. For example a NORM_TRI3 with 3 times the same node id. + */ +void MEDCouplingUMesh::convertDegeneratedCells() throw(INTERP_KERNEL::Exception) +{ + checkFullyDefined(); + if(getMeshDimension()<=1) + throw INTERP_KERNEL::Exception("MEDCouplingUMesh::convertDegeneratedCells works on umeshes with meshdim equals to 2 or 3 !"); + int nbOfCells=getNumberOfCells(); + if(nbOfCells<1) + return ; + int initMeshLgth=getMeshLength(); + int *conn=_nodal_connec->getPointer(); + int *index=_nodal_connec_index->getPointer(); + int posOfCurCell=0; + int newPos=0; + int lgthOfCurCell; + for(int i=0;ireAlloc(newPos); + computeTypes(); +} + /*! * This method checks that all or only polygons (depending 'polyOnly' parameter) 2D cells are correctly oriented relative to 'vec' vector. * The 'vec' vector has to have a non nul norm. @@ -3125,3 +3215,419 @@ void MEDCouplingUMesh::fillInCompact3DMode(int spaceDim, int nbOfNodesInCell, co else throw INTERP_KERNEL::Exception("MEDCouplingUMesh::fillInCompact3DMode : Invalid spaceDim specified : must be 2 or 3 !"); } + + +/*! + * This method takes as input a cell with type 'type' and whose connectivity is defined by (conn,lgth) + * It retrieves the same cell with a potentially different type (in return) whose connectivity is defined by (retConn,retLgth) + * \b WARNING for optimization reason the arrays 'retConn' and 'conn' can overlapped ! + */ +INTERP_KERNEL::NormalizedCellType MEDCouplingUMesh::simplifyDegeneratedCell(INTERP_KERNEL::NormalizedCellType type, const int *conn, int lgth, + int *retConn, int& retLgth) throw(INTERP_KERNEL::Exception) +{ + const INTERP_KERNEL::CellModel& cm=INTERP_KERNEL::CellModel::getCellModel(type); + std::set c(conn,conn+lgth); + c.erase(-1); + bool isObviousNonDegeneratedCell=((int)c.size()==lgth); + if(cm.isQuadratic() || isObviousNonDegeneratedCell) + {//quadratic do nothing for the moment. + retLgth=lgth; + int *tmp=new int[lgth];//no direct std::copy ! overlapping of conn and retConn ! + std::copy(conn,conn+lgth,tmp); + std::copy(tmp,tmp+lgth,retConn); + delete [] tmp; + return type; + } + if(cm.getDimension()==2) + { + int *tmp=new int[lgth]; + tmp[0]=conn[0]; + int newPos=1; + for(int i=1;i faces; + for(unsigned j=0;j nodes(conn,conn+lgth); + nodes.erase(-1); + int nbOfNodes=nodes.size(); + int magicNumber=100*nbOfNodes+nbOfFaces; + switch(magicNumber) + { + case 806: + return tryToUnPolyHex8(conn,nbOfFaces,lgth,retConn,retLgth); + case 506: + return tryToUnPolyPenta6(conn,nbOfFaces,lgth,retConn,retLgth); + case 505: + return tryToUnPolyPyra5(conn,nbOfFaces,lgth,retConn,retLgth); + case 404: + return tryToUnPolyTetra4(conn,nbOfFaces,lgth,retConn,retLgth); + default: + retLgth=lgth; + std::copy(conn,conn+lgth,retConn); + return INTERP_KERNEL::NORM_POLYHED; + } +} + +bool MEDCouplingUMesh::orientOppositeFace(const int *baseFace, int *retConn, const int *sideFace, int lgthBaseFace) +{ + std::vector tmp2; + std::set bases(baseFace,baseFace+lgthBaseFace); + std::set sides(sideFace,sideFace+4); + std::set_intersection(bases.begin(),bases.end(),sides.begin(),sides.end(),std::back_insert_iterator< std::vector >(tmp2)); + if(tmp2.size()!=2) + return false; + std::vector< std::pair > baseEdges(lgthBaseFace); + std::vector< std::pair > oppEdges(lgthBaseFace); + std::vector< std::pair > sideEdges(4); + for(int i=0;i(baseFace[i],baseFace[(i+1)%lgthBaseFace]); + oppEdges[i]=std::pair(retConn[i],retConn[(i+1)%lgthBaseFace]); + } + for(int i=0;i<4;i++) + sideEdges[i]=std::pair(sideFace[i],sideFace[(i+1)%4]); + std::vector< std::pair > tmp; + std::set< std::pair > baseEdgesS(baseEdges.begin(),baseEdges.end()); + std::set< std::pair > sideEdgesS(sideEdges.begin(),sideEdges.end()); + std::set_intersection(baseEdgesS.begin(),baseEdgesS.end(),sideEdgesS.begin(),sideEdgesS.end(),std::back_insert_iterator< std::vector< std::pair > >(tmp)); + if(tmp.empty()) + { + //reverse sideFace + for(int i=0;i<4;i++) + { + std::pair p=sideEdges[i]; + std::pair r(p.second,p.first); + sideEdges[i]=r; + } + //end reverse sideFace + std::set< std::pair > baseEdgesS(baseEdges.begin(),baseEdges.end()); + std::set< std::pair > sideEdgesS(sideEdges.begin(),sideEdges.end()); + std::set_intersection(baseEdgesS.begin(),baseEdgesS.end(),sideEdgesS.begin(),sideEdgesS.end(),std::back_insert_iterator< std::vector< std::pair > >(tmp)); + if(tmp.empty()) + return false; + } + if(tmp.size()!=1) + return false; + bool found=false; + std::pair pInOpp; + for(int i=0;i<4 && !found;i++) + {//finding the pair(edge) in sideFace that do not include any node of tmp[0] edge + found=(tmp[0].first!=sideEdges[i].first && tmp[0].first!=sideEdges[i].second && + tmp[0].second!=sideEdges[i].first && tmp[0].second!=sideEdges[i].second); + if(found) + {//found ! reverse it + pInOpp.first=sideEdges[i].second; + pInOpp.second=sideEdges[i].first; + } + } + if(!found) + return false; + int pos=std::distance(baseEdges.begin(),std::find(baseEdges.begin(),baseEdges.end(),tmp[0])); + std::vector< std::pair >::iterator it=std::find(oppEdges.begin(),oppEdges.end(),pInOpp); + if(it==oppEdges.end())//the opposite edge of side face is not found opposite face ... maybe problem of orientation of polyhedron + return false; + int pos2=std::distance(oppEdges.begin(),it); + int offset=pos-pos2; + if(offset<0) + offset+=lgthBaseFace; + //this is the end copy the result + int *tmp3=new int[lgthBaseFace]; + for(int i=0;i(),(int)INTERP_KERNEL::NORM_QUAD4))==conn+lgth+nbOfFaces) + {//6 faces are QUAD4. + int oppositeFace=-1; + std::set conn1(conn,conn+4); + for(int i=1;i<6 && oppositeFace<0;i++) + { + std::vector tmp; + std::set conn2(conn+5*i,conn+5*i+4); + std::set_intersection(conn1.begin(),conn1.end(),conn2.begin(),conn2.end(),std::back_insert_iterator< std::vector >(tmp)); + if(tmp.empty()) + oppositeFace=i; + } + if(oppositeFace>=1) + {//oppositeFace of face#0 found. + int tmp2[4]; + if(tryToArrangeOppositeFace(conn,lgth,4,conn,conn+5*oppositeFace,6,tmp2)) + { + std::copy(conn,conn+4,retConn); + std::copy(tmp2,tmp2+4,retConn+4); + retLgth=8; + return INTERP_KERNEL::NORM_HEXA8; + } + } + } + retLgth=lgth; + std::copy(conn,conn+lgth,retConn); + return INTERP_KERNEL::NORM_POLYHED; +} + +/*! + * Cell with 'conn' connectivity has been detected as a good candidate. Full check of this. If yes NORM_PENTA6 is returned. + * If fails a POLYHED is returned. + */ +INTERP_KERNEL::NormalizedCellType MEDCouplingUMesh::tryToUnPolyPenta6(const int *conn, int nbOfFaces, int lgth, int *retConn, int& retLgth) +{ + int nbOfTriFace=std::count(conn+lgth,conn+lgth+nbOfFaces,(int)INTERP_KERNEL::NORM_TRI3); + int nbOfQuadFace=std::count(conn+lgth,conn+lgth+nbOfFaces,(int)INTERP_KERNEL::NORM_QUAD4); + if(nbOfTriFace==2 && nbOfQuadFace==3) + { + int tri3_0=std::distance(conn+lgth,std::find(conn+lgth,conn+lgth+nbOfFaces,(int)INTERP_KERNEL::NORM_TRI3)); + int tri3_1=std::distance(conn+lgth,std::find(conn+lgth+tri3_0+1,conn+lgth+nbOfFaces,(int)INTERP_KERNEL::NORM_TRI3)); + const int *tri_0,*tri_1; + const int *w=conn; + for(int i=0;i<5;i++) + { + if(i==tri3_0) + tri_0=w; + if(i==tri3_1) + tri_1=w; + w=std::find(w,conn+lgth,-1); + w++; + } + std::vector tmp; + std::set conn1(tri_0,tri_0+3); + std::set conn2(tri_1,tri_1+3); + std::set_intersection(conn1.begin(),conn1.end(),conn2.begin(),conn2.end(),std::back_insert_iterator< std::vector >(tmp)); + if(tmp.empty()) + { + int tmp2[3]; + if(tryToArrangeOppositeFace(conn,lgth,3,tri_0,tri_1,5,tmp2)) + { + std::copy(conn,conn+4,retConn); + std::copy(tmp2,tmp2+3,retConn+3); + retLgth=6; + return INTERP_KERNEL::NORM_PENTA6; + } + } + } + retLgth=lgth; + std::copy(conn,conn+lgth,retConn); + return INTERP_KERNEL::NORM_POLYHED; +} + +/*! + * Cell with 'conn' connectivity has been detected as a good candidate. Full check of this. If yes NORM_PYRA5 is returned. + * If fails a POLYHED is returned. + */ +INTERP_KERNEL::NormalizedCellType MEDCouplingUMesh::tryToUnPolyPyra5(const int *conn, int nbOfFaces, int lgth, int *retConn, int& retLgth) +{ + int nbOfTriFace=std::count(conn+lgth,conn+lgth+nbOfFaces,(int)INTERP_KERNEL::NORM_TRI3); + int nbOfQuadFace=std::count(conn+lgth,conn+lgth+nbOfFaces,(int)INTERP_KERNEL::NORM_QUAD4); + if(nbOfTriFace==4 && nbOfQuadFace==1) + { + int quad4_pos=std::distance(conn+lgth,std::find(conn+lgth,conn+lgth+nbOfFaces,(int)INTERP_KERNEL::NORM_QUAD4)); + const int *quad4=0; + const int *w=conn; + for(int i=0;i<5 && quad4==0;i++) + { + if(i==quad4_pos) + quad4=w; + w=std::find(w,conn+lgth,-1); + w++; + } + std::set quad4S(quad4,quad4+4); + w=conn; + bool ok=true; + int point=-1; + for(int i=0;i<5 && ok;i++) + if(i!=quad4_pos) + { + std::vector tmp; + std::set conn2(w,w+3); + std::set_intersection(conn2.begin(),conn2.end(),quad4S.begin(),quad4S.end(),std::back_insert_iterator< std::vector >(tmp)); + ok=tmp.size()==2; + tmp.clear(); + std::set_difference(conn2.begin(),conn2.end(),quad4S.begin(),quad4S.end(),std::back_insert_iterator< std::vector >(tmp)); + ok=ok && tmp.size()==1; + if(ok) + { + if(point>=0) + ok=point==tmp[0]; + else + point=tmp[0]; + } + w=std::find(w,conn+lgth,-1); + w++; + } + if(ok && point>=0) + { + std::copy(quad4,quad4+4,retConn); + retConn[4]=point; + retLgth=5; + return INTERP_KERNEL::NORM_PYRA5; + } + } + retLgth=lgth; + std::copy(conn,conn+lgth,retConn); + return INTERP_KERNEL::NORM_POLYHED; +} + +/*! + * Cell with 'conn' connectivity has been detected as a good candidate. Full check of this. If yes NORM_TETRA4 is returned. + * If fails a POLYHED is returned. + */ +INTERP_KERNEL::NormalizedCellType MEDCouplingUMesh::tryToUnPolyTetra4(const int *conn, int nbOfFaces, int lgth, int *retConn, int& retLgth) +{ + if(std::find_if(conn+lgth,conn+lgth+nbOfFaces,std::bind2nd(std::not_equal_to(),(int)INTERP_KERNEL::NORM_TRI3))==conn+lgth+nbOfFaces) + { + std::set tribase(conn,conn+3); + int point=-1; + bool ok=true; + for(int i=1;i<4 && ok;i++) + { + std::vector tmp; + std::set conn2(conn+i*4,conn+4*i+3); + std::set_intersection(conn2.begin(),conn2.end(),tribase.begin(),tribase.end(),std::back_insert_iterator< std::vector >(tmp)); + ok=tmp.size()==2; + tmp.clear(); + std::set_difference(conn2.begin(),conn2.end(),tribase.begin(),tribase.end(),std::back_insert_iterator< std::vector >(tmp)); + ok=ok && tmp.size()==1; + if(ok) + { + if(point>=0) + ok=point==tmp[0]; + else + point=tmp[0]; + } + } + if(ok && point>=0) + { + std::copy(conn,conn+3,retConn); + retConn[3]=point; + retLgth=4; + return INTERP_KERNEL::NORM_TETRA4; + } + } + retLgth=lgth; + std::copy(conn,conn+lgth,retConn); + return INTERP_KERNEL::NORM_POLYHED; +} diff --git a/src/MEDCoupling/MEDCouplingUMesh.hxx b/src/MEDCoupling/MEDCouplingUMesh.hxx index eb989c9a6..7d059f5a1 100644 --- a/src/MEDCoupling/MEDCouplingUMesh.hxx +++ b/src/MEDCoupling/MEDCouplingUMesh.hxx @@ -79,6 +79,7 @@ namespace ParaMEDMEM MEDCOUPLING_EXPORT bool areCellsEqual2(int cell1, int cell2) const; MEDCOUPLING_EXPORT bool areCellsFrom2MeshEqual(const MEDCouplingUMesh *other, int cellId, double prec) const; MEDCOUPLING_EXPORT void convertToPolyTypes(const std::vector& cellIdsToConvert); + MEDCOUPLING_EXPORT void unPolyze(); MEDCOUPLING_EXPORT DataArrayInt *zipCoordsTraducer(); MEDCOUPLING_EXPORT DataArrayInt *zipConnectivityTraducer(int compType); MEDCOUPLING_EXPORT void getReverseNodalConnectivity(DataArrayInt *revNodal, DataArrayInt *revNodalIndx) const; @@ -108,6 +109,7 @@ namespace ParaMEDMEM MEDCOUPLING_EXPORT bool isFullyQuadratic() const; MEDCOUPLING_EXPORT bool isPresenceOfQuadratic() const; MEDCOUPLING_EXPORT void convertQuadraticCellsToLinear() throw(INTERP_KERNEL::Exception); + MEDCOUPLING_EXPORT void convertDegeneratedCells() throw(INTERP_KERNEL::Exception); MEDCOUPLING_EXPORT void are2DCellsNotCorrectlyOriented(const double *vec, bool polyOnly, std::vector& cells) const throw(INTERP_KERNEL::Exception); MEDCOUPLING_EXPORT void orientCorrectly2DCells(const double *vec, bool polyOnly) throw(INTERP_KERNEL::Exception); MEDCOUPLING_EXPORT void arePolyhedronsNotCorrectlyOriented(std::vector& cells) const throw(INTERP_KERNEL::Exception); @@ -155,6 +157,19 @@ namespace ParaMEDMEM double eps, std::vector& elts, std::vector& eltsIndex) const; static void fillInCompact3DMode(int spaceDim, int nbOfNodesInCell, const int *conn, const double *coo, double *zipFrmt) throw(INTERP_KERNEL::Exception); static void appendExtrudedCell(const int *connBg, const int *connEnd, int nbOfNodesPerLev, bool isQuad, std::vector& ret); + static INTERP_KERNEL::NormalizedCellType simplifyDegeneratedCell(INTERP_KERNEL::NormalizedCellType type, const int *conn, int lgth, + int *retConn, int& retLgth) throw(INTERP_KERNEL::Exception); + static int *getFullPolyh3DCell(INTERP_KERNEL::NormalizedCellType type, const int *conn, int lgth, + int& retNbOfFaces, int& retLgth); + static INTERP_KERNEL::NormalizedCellType tryToUnPoly2D(const int *conn, int lgth, int *retConn, int& retLgth); + static INTERP_KERNEL::NormalizedCellType tryToUnPoly3D(const int *conn, int nbOfFaces, int lgth, int *retConn, int& retLgth); + static INTERP_KERNEL::NormalizedCellType tryToUnPolyHex8(const int *conn, int nbOfFaces, int lgth, int *retConn, int& retLgth); + static INTERP_KERNEL::NormalizedCellType tryToUnPolyPenta6(const int *conn, int nbOfFaces, int lgth, int *retConn, int& retLgth); + static INTERP_KERNEL::NormalizedCellType tryToUnPolyPyra5(const int *conn, int nbOfFaces, int lgth, int *retConn, int& retLgth); + static INTERP_KERNEL::NormalizedCellType tryToUnPolyTetra4(const int *conn, int nbOfFaces, int lgth, int *retConn, int& retLgth); + static bool tryToArrangeOppositeFace(const int *conn, int lgth, int lgthBaseFace, const int *baseFace, const int *oppFaceId, int nbOfFaces, int *retConnOfOppFace); + static bool isWellOriented(const int *baseFace, int *retConn, const int *sideFace, int lgthBaseFace); + static bool orientOppositeFace(const int *baseFace, int *retConn, const int *sideFace, int lgthBaseFace); private: //! this iterator stores current position in _nodal_connec array. mutable int _iterator;