X-Git-Url: http://git.salome-platform.org/gitweb/?a=blobdiff_plain;f=src%2FMEDCoupling%2FMEDCouplingUMesh.cxx;h=697b876ab41c2ccf18c9ebb0565217f4e7f06206;hb=fc22b4cd63404700f73e09be7507a11b838e55bc;hp=ce3bb6436b6742dfe93d7b8f38424d1cee3be5fc;hpb=8d510346a2b5a5ec2db2a475eb818f50bff1bb07;p=tools%2Fmedcoupling.git diff --git a/src/MEDCoupling/MEDCouplingUMesh.cxx b/src/MEDCoupling/MEDCouplingUMesh.cxx index ce3bb6436..697b876ab 100644 --- a/src/MEDCoupling/MEDCouplingUMesh.cxx +++ b/src/MEDCoupling/MEDCouplingUMesh.cxx @@ -56,8 +56,8 @@ using namespace MEDCoupling; double MEDCouplingUMesh::EPS_FOR_POLYH_ORIENTATION=1.e-14; /// @cond INTERNAL -const INTERP_KERNEL::NormalizedCellType MEDCouplingUMesh::MEDMEM_ORDER[N_MEDMEM_ORDER] = { INTERP_KERNEL::NORM_POINT1, INTERP_KERNEL::NORM_SEG2, INTERP_KERNEL::NORM_SEG3, INTERP_KERNEL::NORM_SEG4, INTERP_KERNEL::NORM_POLYL, INTERP_KERNEL::NORM_TRI3, INTERP_KERNEL::NORM_QUAD4, INTERP_KERNEL::NORM_TRI6, INTERP_KERNEL::NORM_TRI7, INTERP_KERNEL::NORM_QUAD8, INTERP_KERNEL::NORM_QUAD9, INTERP_KERNEL::NORM_POLYGON, INTERP_KERNEL::NORM_QPOLYG, INTERP_KERNEL::NORM_TETRA4, INTERP_KERNEL::NORM_PYRA5, INTERP_KERNEL::NORM_PENTA6, INTERP_KERNEL::NORM_HEXA8, INTERP_KERNEL::NORM_HEXGP12, INTERP_KERNEL::NORM_TETRA10, INTERP_KERNEL::NORM_PYRA13, INTERP_KERNEL::NORM_PENTA15, INTERP_KERNEL::NORM_HEXA20, INTERP_KERNEL::NORM_HEXA27, INTERP_KERNEL::NORM_POLYHED }; -const int MEDCouplingUMesh::MEDCOUPLING2VTKTYPETRADUCER[INTERP_KERNEL::NORM_MAXTYPE+1]={1,3,21,5,9,7,22,34,23,28,-1,-1,-1,-1,10,14,13,-1,12,-1,24,-1,16,27,-1,26,-1,29,-1,-1,25,42,36,4}; +const INTERP_KERNEL::NormalizedCellType MEDCouplingUMesh::MEDMEM_ORDER[N_MEDMEM_ORDER] = { INTERP_KERNEL::NORM_POINT1, INTERP_KERNEL::NORM_SEG2, INTERP_KERNEL::NORM_SEG3, INTERP_KERNEL::NORM_SEG4, INTERP_KERNEL::NORM_POLYL, INTERP_KERNEL::NORM_TRI3, INTERP_KERNEL::NORM_QUAD4, INTERP_KERNEL::NORM_TRI6, INTERP_KERNEL::NORM_TRI7, INTERP_KERNEL::NORM_QUAD8, INTERP_KERNEL::NORM_QUAD9, INTERP_KERNEL::NORM_POLYGON, INTERP_KERNEL::NORM_QPOLYG, INTERP_KERNEL::NORM_TETRA4, INTERP_KERNEL::NORM_PYRA5, INTERP_KERNEL::NORM_PENTA6, INTERP_KERNEL::NORM_HEXA8, INTERP_KERNEL::NORM_HEXGP12, INTERP_KERNEL::NORM_TETRA10, INTERP_KERNEL::NORM_PYRA13, INTERP_KERNEL::NORM_PENTA15, INTERP_KERNEL::NORM_PENTA18, INTERP_KERNEL::NORM_HEXA20, INTERP_KERNEL::NORM_HEXA27, INTERP_KERNEL::NORM_POLYHED }; +const int MEDCouplingUMesh::MEDCOUPLING2VTKTYPETRADUCER[INTERP_KERNEL::NORM_MAXTYPE+1]={1,3,21,5,9,7,22,34,23,28,-1,-1,-1,-1,10,14,13,-1,12,-1,24,-1,16,27,-1,26,-1,29,32,-1,25,42,36,4}; /// @endcond MEDCouplingUMesh *MEDCouplingUMesh::New() @@ -1984,19 +1984,19 @@ void MEDCouplingUMesh::setPartOfMySelf(const int *cellIdsBg, const int *cellIdsE oss << ", whereas other mesh dimension is set equal to " << otherOnSameCoordsThanThis.getMeshDimension() << " !"; throw INTERP_KERNEL::Exception(oss.str()); } - int nbOfCellsToModify=(int)std::distance(cellIdsBg,cellIdsEnd); + std::size_t nbOfCellsToModify(std::distance(cellIdsBg,cellIdsEnd)); if(nbOfCellsToModify!=otherOnSameCoordsThanThis.getNumberOfCells()) { std::ostringstream oss; oss << "MEDCouplingUMesh::setPartOfMySelf : cells ids length (" << nbOfCellsToModify << ") do not match the number of cells of other mesh (" << otherOnSameCoordsThanThis.getNumberOfCells() << ") !"; throw INTERP_KERNEL::Exception(oss.str()); } - int nbOfCells=getNumberOfCells(); - bool easyAssign=true; - const int *connI=_nodal_connec_index->getConstPointer(); - const int *connIOther=otherOnSameCoordsThanThis._nodal_connec_index->getConstPointer(); + std::size_t nbOfCells(getNumberOfCells()); + bool easyAssign(true); + const int *connI(_nodal_connec_index->begin()); + const int *connIOther=otherOnSameCoordsThanThis._nodal_connec_index->begin(); for(const int *it=cellIdsBg;it!=cellIdsEnd && easyAssign;it++,connIOther++) { - if(*it>=0 && *it=0 && *it<(int)nbOfCells) { easyAssign=(connIOther[1]-connIOther[0])==(connI[*it+1]-connI[*it]); } @@ -2034,7 +2034,7 @@ void MEDCouplingUMesh::setPartOfMySelfSlice(int start, int end, int step, const throw INTERP_KERNEL::Exception(oss.str()); } int nbOfCellsToModify=DataArray::GetNumberOfItemGivenBESRelative(start,end,step,"MEDCouplingUMesh::setPartOfMySelfSlice : "); - if(nbOfCellsToModify!=otherOnSameCoordsThanThis.getNumberOfCells()) + if(nbOfCellsToModify!=(int)otherOnSameCoordsThanThis.getNumberOfCells()) { std::ostringstream oss; oss << "MEDCouplingUMesh::setPartOfMySelfSlice : cells ids length (" << nbOfCellsToModify << ") do not match the number of cells of other mesh (" << otherOnSameCoordsThanThis.getNumberOfCells() << ") !"; throw INTERP_KERNEL::Exception(oss.str()); @@ -2363,7 +2363,7 @@ void MEDCouplingUMesh::findNodesToDuplicate(const MEDCouplingUMesh& otherDimM1On DAInt neighIInit00(tmp11); // Neighbor information of the mesh WITH the crack (some neighbors are removed): DataArrayInt *idsTmp=0; - bool b=m01->areCellsIncludedIn(&otherDimM1OnSameCoords,2,idsTmp); + m01->areCellsIncludedIn(&otherDimM1OnSameCoords,2,idsTmp); DAInt ids(idsTmp); // In the neighbor information remove the connection between high dimension cells and its low level constituents which are part // of the frontier given in parameter (i.e. the cells of low dimension from the group delimiting the crack): @@ -2390,7 +2390,7 @@ void MEDCouplingUMesh::findNodesToDuplicate(const MEDCouplingUMesh& otherDimM1On // Connex zone without the crack (to compute the next seed really) int dnu; DAInt connexCheck = MEDCouplingUMesh::ComputeSpreadZoneGraduallyFromSeed(&seed, &seed+1, neighInit00,neighIInit00, -1, dnu); - int cnt = 0; + std::size_t cnt(0); for (int * ptr = connexCheck->getPointer(); cnt < connexCheck->getNumberOfTuples(); ptr++, cnt++) hitCells->setIJ(*ptr,0,1); // Connex zone WITH the crack (to identify cells lying on either part of the crack) @@ -2778,7 +2778,7 @@ DataArrayInt *MEDCouplingUMesh::getCellsInBoundingBox(const INTERP_KERNEL::Direc INTERP_KERNEL::NormalizedCellType MEDCouplingUMesh::getTypeOfCell(std::size_t cellId) const { const int *ptI(_nodal_connec_index->begin()),*pt(_nodal_connec->begin()); - if(cellId>=0 && cellId<_nodal_connec_index->getNbOfElems()-1) + if(cellId<_nodal_connec_index->getNbOfElems()-1) return (INTERP_KERNEL::NormalizedCellType) pt[ptI[cellId]]; else { @@ -2839,10 +2839,9 @@ std::size_t MEDCouplingUMesh::getNumberOfCellsWithType(INTERP_KERNEL::Normalized * cleared before the appending. * \throw If \a cellId is invalid. Valid range is [0, \a this->getNumberOfCells() ). */ -void MEDCouplingUMesh::getNodeIdsOfCell(int cellId, std::vector& conn) const +void MEDCouplingUMesh::getNodeIdsOfCell(std::size_t cellId, std::vector& conn) const { - const int *ptI=_nodal_connec_index->getConstPointer(); - const int *pt=_nodal_connec->getConstPointer(); + const int *ptI(_nodal_connec_index->begin()),*pt(_nodal_connec->begin()); for(const int *w=pt+ptI[cellId]+1;w!=pt+ptI[cellId+1];w++) if(*w>=0) conn.push_back(*w); @@ -3041,14 +3040,14 @@ void MEDCouplingUMesh::setConnectivity(DataArrayInt *conn, DataArrayInt *connInd * Copy constructor. If 'deepCopy' is false \a this is a shallow copy of other. * If 'deeCpy' is true all arrays (coordinates and connectivities) are deeply copied. */ -MEDCouplingUMesh::MEDCouplingUMesh(const MEDCouplingUMesh& other, bool deepCopy):MEDCouplingPointSet(other,deepCopy),_mesh_dim(other._mesh_dim), +MEDCouplingUMesh::MEDCouplingUMesh(const MEDCouplingUMesh& other, bool deepCpy):MEDCouplingPointSet(other,deepCpy),_mesh_dim(other._mesh_dim), _nodal_connec(0),_nodal_connec_index(0), _types(other._types) { if(other._nodal_connec) - _nodal_connec=other._nodal_connec->performCopyOrIncrRef(deepCopy); + _nodal_connec=other._nodal_connec->performCopyOrIncrRef(deepCpy); if(other._nodal_connec_index) - _nodal_connec_index=other._nodal_connec_index->performCopyOrIncrRef(deepCopy); + _nodal_connec_index=other._nodal_connec_index->performCopyOrIncrRef(deepCpy); } MEDCouplingUMesh::~MEDCouplingUMesh() @@ -3879,7 +3878,7 @@ void MEDCouplingUMesh::project1D(const double *pt, const double *v, double eps, MCAuto f=buildDirectionVectorField(); const double *fPtr=f->getArray()->getConstPointer(); double tmp[3]; - for(int i=0;igetNumberOfComponents()!=spaceDim) + if((int)pts->getNumberOfComponents()!=spaceDim) { std::ostringstream oss; oss << "MEDCouplingUMesh::distanceToPoints : input pts DataArrayDouble has " << pts->getNumberOfComponents() << " components whereas it should be equal to " << spaceDim << " (mesh spaceDimension) !"; throw INTERP_KERNEL::Exception(oss.str()); @@ -5418,6 +5417,8 @@ DataArrayDouble *MEDCouplingUMesh::getBoundingBoxForBBTreeFast() const DataArrayDouble *MEDCouplingUMesh::getBoundingBoxForBBTree2DQuadratic(double arcDetEps) const { checkFullyDefined(); + INTERP_KERNEL::QuadraticPlanarArcDetectionPrecision arcPrec(arcDetEps); + int spaceDim(getSpaceDimension()),mDim(getMeshDimension()),nbOfCells(getNumberOfCells()); if(spaceDim!=2 || mDim!=2) throw INTERP_KERNEL::Exception("MEDCouplingUMesh::getBoundingBoxForBBTree2DQuadratic : This method should be applied on mesh with mesh dimension equal to 2 and space dimension also equal to 2!"); @@ -5429,7 +5430,6 @@ DataArrayDouble *MEDCouplingUMesh::getBoundingBoxForBBTree2DQuadratic(double arc { const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel((INTERP_KERNEL::NormalizedCellType)conn[*connI])); int sz(connI[1]-connI[0]-1); - INTERP_KERNEL::QUADRATIC_PLANAR::_arc_detection_precision=arcDetEps; std::vector nodes(sz); INTERP_KERNEL::QuadraticPolygon *pol(0); for(int j=0;j ret(DataArrayDouble::New()); ret->alloc(nbOfCells,2*spaceDim); double *bbox(ret->getPointer()); const double *coords(_coords->begin()); @@ -5474,7 +5475,6 @@ DataArrayDouble *MEDCouplingUMesh::getBoundingBoxForBBTree1DQuadratic(double arc { const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel((INTERP_KERNEL::NormalizedCellType)conn[*connI])); int sz(connI[1]-connI[0]-1); - INTERP_KERNEL::QUADRATIC_PLANAR::_arc_detection_precision=arcDetEps; std::vector nodes(sz); INTERP_KERNEL::Edge *edge(0); for(int j=0;j=3) + if(nodalI[1]-nodalI[0]>=4) { + double aa[3]={coor[nodal[nodalI[0]+1+1]*3+0]-coor[nodal[nodalI[0]+1+0]*3+0], + coor[nodal[nodalI[0]+1+1]*3+1]-coor[nodal[nodalI[0]+1+0]*3+1], + coor[nodal[nodalI[0]+1+1]*3+2]-coor[nodal[nodalI[0]+1+0]*3+2]} + ,bb[3]={coor[nodal[nodalI[0]+1+2]*3+0]-coor[nodal[nodalI[0]+1+0]*3+0], + coor[nodal[nodalI[0]+1+2]*3+1]-coor[nodal[nodalI[0]+1+0]*3+1], + coor[nodal[nodalI[0]+1+2]*3+2]-coor[nodal[nodalI[0]+1+0]*3+2]}; + double cc[3]={aa[1]*bb[2]-aa[2]*bb[1],aa[2]*bb[0]-aa[0]*bb[2],aa[0]*bb[1]-aa[1]*bb[0]}; for(int j=0;j<3;j++) { int nodeId(nodal[nodalI[0]+1+j]); @@ -6489,14 +6496,34 @@ DataArrayDouble *MEDCouplingUMesh::computePlaneEquationOf3DFaces() const throw INTERP_KERNEL::Exception(oss.str()); } } + if(sqrt(cc[0]*cc[0]+cc[1]*cc[1]+cc[2]*cc[2])>1e-7) + { + INTERP_KERNEL::inverseMatrix(matrix,4,matrix2); + retPtr[0]=matrix2[3]; retPtr[1]=matrix2[7]; retPtr[2]=matrix2[11]; retPtr[3]=matrix2[15]; + } + else + { + if(nodalI[1]-nodalI[0]==4) + { + std::ostringstream oss; oss << "MEDCouplingUMesh::computePlaneEquationOf3DFaces : cell" << i << " : Presence of The 3 colinear points !"; + throw INTERP_KERNEL::Exception(oss.str()); + } + // + double dd[3]={0.,0.,0.}; + for(int offset=nodalI[0]+1;offset()); + int nbOfNodesInCell(nodalI[1]-nodalI[0]-1); + std::transform(dd,dd+3,dd,std::bind2nd(std::multiplies(),1./(double)nbOfNodesInCell)); + std::copy(dd,dd+3,matrix+4*2); + INTERP_KERNEL::inverseMatrix(matrix,4,matrix2); + retPtr[0]=matrix2[3]; retPtr[1]=matrix2[7]; retPtr[2]=matrix2[11]; retPtr[3]=matrix2[15]; + } } else { std::ostringstream oss; oss << "MEDCouplingUMesh::computePlaneEquationOf3DFaces : invalid 2D cell #" << i << " ! Must be constitued by more than 3 nodes !"; throw INTERP_KERNEL::Exception(oss.str()); } - INTERP_KERNEL::inverseMatrix(matrix,4,matrix2); - retPtr[0]=matrix2[3]; retPtr[1]=matrix2[7]; retPtr[2]=matrix2[11]; retPtr[3]=matrix2[15]; } return ret.retn(); } @@ -7700,7 +7727,7 @@ bool MEDCouplingUMesh::RemoveIdsFromIndexedArrays(const int *idsToRemoveBg, cons previousArrI=*arrIPtr; *arrIPtr=(int)arrOut.size(); } - if(arr->getNumberOfTuples()==(int)arrOut.size()) + if(arr->getNumberOfTuples()==arrOut.size()) return false; arr->alloc((int)arrOut.size(),1); std::copy(arrOut.begin(),arrOut.end(),arr->getPointer());