X-Git-Url: http://git.salome-platform.org/gitweb/?a=blobdiff_plain;f=src%2FMEDCoupling%2FMEDCouplingUMesh.cxx;h=7a4aae408c4ec5742cdfd77c61269c1eab2ead3c;hb=659f8c67d0348350e12fde38fe8c4de1ff95dffe;hp=88cbe90bc6d630a1594c1a5decaf4a389ab5d56b;hpb=5789a86f2ce485bfafd5740f1ecfa140f9a0e481;p=tools%2Fmedcoupling.git diff --git a/src/MEDCoupling/MEDCouplingUMesh.cxx b/src/MEDCoupling/MEDCouplingUMesh.cxx index 88cbe90bc..7a4aae408 100644 --- a/src/MEDCoupling/MEDCouplingUMesh.cxx +++ b/src/MEDCoupling/MEDCouplingUMesh.cxx @@ -1,9 +1,9 @@ -// Copyright (C) 2007-2013 CEA/DEN, EDF R&D +// Copyright (C) 2007-2014 CEA/DEN, EDF R&D // // This library is free software; you can redistribute it and/or // modify it under the terms of the GNU Lesser General Public // License as published by the Free Software Foundation; either -// version 2.1 of the License. +// version 2.1 of the License, or (at your option) any later version. // // This library is distributed in the hope that it will be useful, // but WITHOUT ANY WARRANTY; without even the implied warranty of @@ -30,10 +30,10 @@ #include "BBTreeDst.txx" #include "SplitterTetra.hxx" #include "DirectedBoundingBox.hxx" +#include "InterpKernelMatrixTools.hxx" #include "InterpKernelMeshQuality.hxx" #include "InterpKernelCellSimplify.hxx" #include "InterpKernelGeo2DEdgeArcCircle.hxx" -#include "MEDCouplingAutoRefCountObjectPtr.hxx" #include "InterpKernelAutoPtr.hxx" #include "InterpKernelGeo2DNode.hxx" #include "InterpKernelGeo2DEdgeLin.hxx" @@ -58,7 +58,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::New() return new MEDCouplingUMesh; } -MEDCouplingUMesh *MEDCouplingUMesh::New(const char *meshName, int meshDim) +MEDCouplingUMesh *MEDCouplingUMesh::New(const std::string& meshName, int meshDim) { MEDCouplingUMesh *ret=new MEDCouplingUMesh; ret->setName(meshName); @@ -96,7 +96,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::clone(bool recDeepCpy) const * \return MEDCouplingUMesh * - A new object instance holding the copy of \a this (deep for connectivity, shallow for coordiantes) * \sa MEDCouplingUMesh::deepCpy */ -MEDCouplingPointSet *MEDCouplingUMesh::deepCpyConnectivityOnly() const throw(INTERP_KERNEL::Exception) +MEDCouplingPointSet *MEDCouplingUMesh::deepCpyConnectivityOnly() const { checkConnectivityFullyDefined(); MEDCouplingAutoRefCountObjectPtr ret=clone(false); @@ -105,7 +105,7 @@ MEDCouplingPointSet *MEDCouplingUMesh::deepCpyConnectivityOnly() const throw(INT return ret.retn(); } -void MEDCouplingUMesh::shallowCopyConnectivityFrom(const MEDCouplingPointSet *other) throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::shallowCopyConnectivityFrom(const MEDCouplingPointSet *other) { if(!other) throw INTERP_KERNEL::Exception("MEDCouplingUMesh::shallowCopyConnectivityFrom : input pointer is null !"); @@ -116,14 +116,20 @@ void MEDCouplingUMesh::shallowCopyConnectivityFrom(const MEDCouplingPointSet *ot setConnectivity(otherC2->getNodalConnectivity(),otherC2->getNodalConnectivityIndex(),true); } -std::size_t MEDCouplingUMesh::getHeapMemorySize() const +std::size_t MEDCouplingUMesh::getHeapMemorySizeWithoutChildren() const { - std::size_t ret=0; + std::size_t ret(MEDCouplingPointSet::getHeapMemorySizeWithoutChildren()); + return ret; +} + +std::vector MEDCouplingUMesh::getDirectChildren() const +{ + std::vector ret(MEDCouplingPointSet::getDirectChildren()); if(_nodal_connec) - ret+=_nodal_connec->getHeapMemorySize(); + ret.push_back(_nodal_connec); if(_nodal_connec_index) - ret+=_nodal_connec_index->getHeapMemorySize(); - return MEDCouplingPointSet::getHeapMemorySize()+ret; + ret.push_back(_nodal_connec_index); + return ret; } void MEDCouplingUMesh::updateTime() const @@ -156,7 +162,7 @@ MEDCouplingUMesh::MEDCouplingUMesh():_mesh_dim(-2),_nodal_connec(0),_nodal_conne * \throw If the connectivity index data array has more than one component. * \throw If the connectivity index data array has a named component. */ -void MEDCouplingUMesh::checkCoherency() const throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::checkCoherency() const { if(_mesh_dim<-1) throw INTERP_KERNEL::Exception("No mesh dimension specified !"); @@ -209,7 +215,7 @@ void MEDCouplingUMesh::checkCoherency() const throw(INTERP_KERNEL::Exception) * \throw If number of nodes defining an element does not correspond to the type of element. * \throw If the nodal connectivity includes an invalid node id. */ -void MEDCouplingUMesh::checkCoherency1(double eps) const throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::checkCoherency1(double eps) const { checkCoherency(); if(_mesh_dim==-1) @@ -281,7 +287,7 @@ void MEDCouplingUMesh::checkCoherency1(double eps) const throw(INTERP_KERNEL::Ex * \throw If number of nodes defining an element does not correspond to the type of element. * \throw If the nodal connectivity includes an invalid node id. */ -void MEDCouplingUMesh::checkCoherency2(double eps) const throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::checkCoherency2(double eps) const { checkCoherency1(eps); } @@ -342,7 +348,7 @@ void MEDCouplingUMesh::allocateCells(int nbOfCells) * \ref medcouplingcppexamplesUmeshStdBuild1 "Here is a C++ example".
* \ref medcouplingpyexamplesUmeshStdBuild1 "Here is a Python example". */ -void MEDCouplingUMesh::insertNextCell(INTERP_KERNEL::NormalizedCellType type, int size, const int *nodalConnOfCell) throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::insertNextCell(INTERP_KERNEL::NormalizedCellType type, int size, const int *nodalConnOfCell) { const INTERP_KERNEL::CellModel& cm=INTERP_KERNEL::CellModel::GetCellModel(type); if(_nodal_connec_index==0) @@ -402,7 +408,7 @@ MEDCouplingUMeshCellIterator *MEDCouplingUMesh::cellIterator() * In this case MEDCouplingUMesh::sortCellsInMEDFileFrmt or MEDCouplingUMesh::rearrange2ConsecutiveCellTypes methods for example can be called before invoking this method. * Useful for python users. */ -MEDCouplingUMeshCellByTypeEntry *MEDCouplingUMesh::cellsByType() throw(INTERP_KERNEL::Exception) +MEDCouplingUMeshCellByTypeEntry *MEDCouplingUMesh::cellsByType() { if(!checkConsecutiveCellTypes()) throw INTERP_KERNEL::Exception("MEDCouplingUMesh::cellsByType : this mesh is not sorted by type !"); @@ -413,17 +419,44 @@ MEDCouplingUMeshCellByTypeEntry *MEDCouplingUMesh::cellsByType() throw(INTERP_KE * Returns a set of all cell types available in \a this mesh. * \return std::set - the set of cell types. * \warning this method does not throw any exception even if \a this is not defined. + * \sa MEDCouplingUMesh::getAllGeoTypesSorted */ std::set MEDCouplingUMesh::getAllGeoTypes() const { return _types; } +/*! + * This method returns the sorted list of geometric types in \a this. + * Sorted means in the same order than the cells in \a this. A single entry in return vector means the maximal chunk of consecutive cells in \a this + * having the same geometric type. So a same geometric type can appear more than once if the cells are not sorted per geometric type. + * + * \throw if connectivity in \a this is not correctly defined. + * + * \sa MEDCouplingMesh::getAllGeoTypes + */ +std::vector MEDCouplingUMesh::getAllGeoTypesSorted() const +{ + std::vector ret; + checkConnectivityFullyDefined(); + int nbOfCells(getNumberOfCells()); + if(nbOfCells==0) + return ret; + if(getMeshLength()<1) + throw INTERP_KERNEL::Exception("MEDCouplingUMesh::getAllGeoTypesSorted : the connectivity in this seems invalid !"); + const int *c(_nodal_connec->begin()),*ci(_nodal_connec_index->begin()); + ret.push_back((INTERP_KERNEL::NormalizedCellType)c[*ci++]); + for(int i=1;i(other); @@ -553,7 +586,7 @@ void MEDCouplingUMesh::checkFastEquivalWith(const MEDCouplingMesh *other, double * \ref cpp_mcumesh_getReverseNodalConnectivity "Here is a C++ example".
* \ref py_mcumesh_getReverseNodalConnectivity "Here is a Python example". */ -void MEDCouplingUMesh::getReverseNodalConnectivity(DataArrayInt *revNodal, DataArrayInt *revNodalIndx) const throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::getReverseNodalConnectivity(DataArrayInt *revNodal, DataArrayInt *revNodalIndx) const { checkFullyDefined(); int nbOfNodes=getNumberOfNodes(); @@ -686,7 +719,7 @@ private: * \ref py_mcumesh_buildDescendingConnectivity "Here is a Python example". * \sa buildDescendingConnectivity2() */ -MEDCouplingUMesh *MEDCouplingUMesh::buildDescendingConnectivity(DataArrayInt *desc, DataArrayInt *descIndx, DataArrayInt *revDesc, DataArrayInt *revDescIndx) const throw(INTERP_KERNEL::Exception) +MEDCouplingUMesh *MEDCouplingUMesh::buildDescendingConnectivity(DataArrayInt *desc, DataArrayInt *descIndx, DataArrayInt *revDesc, DataArrayInt *revDescIndx) const { return buildDescendingConnectivityGen(desc,descIndx,revDesc,revDescIndx,MEDCouplingFastNbrer); } @@ -698,7 +731,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::buildDescendingConnectivity(DataArrayInt *de * This method returns 4 arrays and a mesh as MEDCouplingUMesh::buildDescendingConnectivity does. * \sa MEDCouplingUMesh::buildDescendingConnectivity */ -MEDCouplingUMesh *MEDCouplingUMesh::explode3DMeshTo1D(DataArrayInt *desc, DataArrayInt *descIndx, DataArrayInt *revDesc, DataArrayInt *revDescIndx) const throw(INTERP_KERNEL::Exception) +MEDCouplingUMesh *MEDCouplingUMesh::explode3DMeshTo1D(DataArrayInt *desc, DataArrayInt *descIndx, DataArrayInt *revDesc, DataArrayInt *revDescIndx) const { checkFullyDefined(); if(getMeshDimension()!=3) @@ -755,7 +788,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::explode3DMeshTo1D(DataArrayInt *desc, DataAr * \ref py_mcumesh_buildDescendingConnectivity2 "Here is a Python example". * \sa buildDescendingConnectivity() */ -MEDCouplingUMesh *MEDCouplingUMesh::buildDescendingConnectivity2(DataArrayInt *desc, DataArrayInt *descIndx, DataArrayInt *revDesc, DataArrayInt *revDescIndx) const throw(INTERP_KERNEL::Exception) +MEDCouplingUMesh *MEDCouplingUMesh::buildDescendingConnectivity2(DataArrayInt *desc, DataArrayInt *descIndx, DataArrayInt *revDesc, DataArrayInt *revDescIndx) const { return buildDescendingConnectivityGen(desc,descIndx,revDesc,revDescIndx,MEDCouplingOrientationSensitiveNbrer); } @@ -770,7 +803,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::buildDescendingConnectivity2(DataArrayInt *d * parameter allows to select the right part in this array. The number of tuples is equal to the last values in \b neighborsIndx. * \param [out] neighborsIndx is an array of size this->getNumberOfCells()+1 newly allocated and should be dealt by the caller. This arrays allow to use the first output parameter \b neighbors. */ -void MEDCouplingUMesh::computeNeighborsOfCells(DataArrayInt *&neighbors, DataArrayInt *&neighborsIndx) const throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::computeNeighborsOfCells(DataArrayInt *&neighbors, DataArrayInt *&neighborsIndx) const { MEDCouplingAutoRefCountObjectPtr desc=DataArrayInt::New(); MEDCouplingAutoRefCountObjectPtr descIndx=DataArrayInt::New(); @@ -834,7 +867,7 @@ void MEDCouplingUMesh::ComputeNeighborsOfCellsAdv(const DataArrayInt *desc, cons * For speed reasons no check of this will be done. */ template -MEDCouplingUMesh *MEDCouplingUMesh::buildDescendingConnectivityGen(DataArrayInt *desc, DataArrayInt *descIndx, DataArrayInt *revDesc, DataArrayInt *revDescIndx, DimM1DescNbrer nbrer) const throw(INTERP_KERNEL::Exception) +MEDCouplingUMesh *MEDCouplingUMesh::buildDescendingConnectivityGen(DataArrayInt *desc, DataArrayInt *descIndx, DataArrayInt *revDesc, DataArrayInt *revDescIndx, DimM1DescNbrer nbrer) const { if(!desc || !descIndx || !revDesc || !revDescIndx) throw INTERP_KERNEL::Exception("MEDCouplingUMesh::buildDescendingConnectivityGen : present of a null pointer in input !"); @@ -846,7 +879,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::buildDescendingConnectivityGen(DataArrayInt const int *conn=_nodal_connec->getConstPointer(); const int *connIndex=_nodal_connec_index->getConstPointer(); std::string name="Mesh constituent of "; name+=getName(); - MEDCouplingAutoRefCountObjectPtr ret=MEDCouplingUMesh::New(name.c_str(),getMeshDimension()-SonsGenerator::DELTA); + MEDCouplingAutoRefCountObjectPtr ret=MEDCouplingUMesh::New(name,getMeshDimension()-SonsGenerator::DELTA); ret->setCoords(getCoords()); ret->allocateCells(2*nbOfCells); descIndx->alloc(nbOfCells+1,1); @@ -1096,7 +1129,7 @@ void MEDCouplingUMesh::convertAllToPoly() * \ref cpp_mcumesh_arePolyhedronsNotCorrectlyOriented "Here is a C++ example".
* \ref py_mcumesh_arePolyhedronsNotCorrectlyOriented "Here is a Python example". */ -void MEDCouplingUMesh::convertExtrudedPolyhedra() throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::convertExtrudedPolyhedra() { checkFullyDefined(); if(getMeshDimension()!=3 || getSpaceDimension()!=3) @@ -1250,7 +1283,7 @@ bool MEDCouplingUMesh::unPolyze() * \param [in] eps is a relative precision that allows to establish if some 3D plane are coplanar or not. This epsilon is used to recenter around origin to have maximal * precision. */ -void MEDCouplingUMesh::simplifyPolyhedra(double eps) throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::simplifyPolyhedra(double eps) { checkFullyDefined(); if(getMeshDimension()!=3 || getSpaceDimension()!=3) @@ -1289,7 +1322,7 @@ void MEDCouplingUMesh::simplifyPolyhedra(double eps) throw(INTERP_KERNEL::Except * \return a newly allocated DataArrayInt sorted ascendingly of fetched node ids. * \sa MEDCouplingUMesh::getNodeIdsInUse */ -DataArrayInt *MEDCouplingUMesh::computeFetchedNodeIds() const throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::computeFetchedNodeIds() const { checkConnectivityFullyDefined(); int nbOfCells=getNumberOfCells(); @@ -1319,7 +1352,7 @@ DataArrayInt *MEDCouplingUMesh::computeFetchedNodeIds() const throw(INTERP_KERNE * \param [in,out] nodeIdsInUse an array of size typically equal to nbOfNodes. * \sa MEDCouplingUMesh::getNodeIdsInUse */ -void MEDCouplingUMesh::computeNodeIdsAlg(std::vector& nodeIdsInUse) const throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::computeNodeIdsAlg(std::vector& nodeIdsInUse) const { int nbOfNodes=(int)nodeIdsInUse.size(); int nbOfCells=getNumberOfCells(); @@ -1356,7 +1389,7 @@ void MEDCouplingUMesh::computeNodeIdsAlg(std::vector& nodeIdsInUse) const * \ref py_mcumesh_getNodeIdsInUse "Here is a Python example". * \sa computeNodeIdsAlg() */ -DataArrayInt *MEDCouplingUMesh::getNodeIdsInUse(int& nbrOfNodesInUse) const throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::getNodeIdsInUse(int& nbrOfNodesInUse) const { nbrOfNodesInUse=-1; int nbOfNodes=getNumberOfNodes(); @@ -1393,7 +1426,7 @@ DataArrayInt *MEDCouplingUMesh::getNodeIdsInUse(int& nbrOfNodesInUse) const thro * \return a newly allocated array * \sa MEDCouplingUMesh::computeEffectiveNbOfNodesPerCell */ -DataArrayInt *MEDCouplingUMesh::computeNbOfNodesPerCell() const throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::computeNbOfNodesPerCell() const { checkConnectivityFullyDefined(); int nbOfCells=getNumberOfCells(); @@ -1419,7 +1452,7 @@ DataArrayInt *MEDCouplingUMesh::computeNbOfNodesPerCell() const throw(INTERP_KER * \return DataArrayInt * - new object to be deallocated by the caller. * \sa MEDCouplingUMesh::computeNbOfNodesPerCell */ -DataArrayInt *MEDCouplingUMesh::computeEffectiveNbOfNodesPerCell() const throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::computeEffectiveNbOfNodesPerCell() const { checkConnectivityFullyDefined(); int nbOfCells=getNumberOfCells(); @@ -1448,7 +1481,7 @@ DataArrayInt *MEDCouplingUMesh::computeEffectiveNbOfNodesPerCell() const throw(I * * \return a newly allocated array */ -DataArrayInt *MEDCouplingUMesh::computeNbOfFacesPerCell() const throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::computeNbOfFacesPerCell() const { checkConnectivityFullyDefined(); int nbOfCells=getNumberOfCells(); @@ -1479,7 +1512,7 @@ DataArrayInt *MEDCouplingUMesh::computeNbOfFacesPerCell() const throw(INTERP_KER * \ref cpp_mcumesh_zipCoordsTraducer "Here is a C++ example".
* \ref py_mcumesh_zipCoordsTraducer "Here is a Python example". */ -DataArrayInt *MEDCouplingUMesh::zipCoordsTraducer() throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::zipCoordsTraducer() { return MEDCouplingPointSet::zipCoordsTraducer(); } @@ -1692,7 +1725,7 @@ bool MEDCouplingUMesh::AreCellsEqualInPool(const std::vector& candidates, i * \return the correspondance array old to new in a newly allocated array. * */ -void MEDCouplingUMesh::findCommonCells(int compType, int startCellId, DataArrayInt *& commonCellsArr, DataArrayInt *& commonCellsIArr) const throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::findCommonCells(int compType, int startCellId, DataArrayInt *& commonCellsArr, DataArrayInt *& commonCellsIArr) const { MEDCouplingAutoRefCountObjectPtr revNodal=DataArrayInt::New(),revNodalI=DataArrayInt::New(); getReverseNodalConnectivity(revNodal,revNodalI); @@ -1802,7 +1835,7 @@ void MEDCouplingUMesh::FindCommonCellsAlg(int compType, int startCellId, const D * \sa checkDeepEquivalOnSameNodesWith() * \sa checkGeoEquivalWith() */ -bool MEDCouplingUMesh::areCellsIncludedIn(const MEDCouplingUMesh *other, int compType, DataArrayInt *& arr) const throw(INTERP_KERNEL::Exception) +bool MEDCouplingUMesh::areCellsIncludedIn(const MEDCouplingUMesh *other, int compType, DataArrayInt *& arr) const { MEDCouplingAutoRefCountObjectPtr mesh=MergeUMeshesOnSameCoords(this,other); int nbOfCells=getNumberOfCells(); @@ -1816,7 +1849,7 @@ bool MEDCouplingUMesh::areCellsIncludedIn(const MEDCouplingUMesh *other, int com } MEDCouplingAutoRefCountObjectPtr o2n=mesh->zipConnectivityTraducer(compType,nbOfCells); arr=o2n->substr(nbOfCells); - arr->setName(other->getName().c_str()); + arr->setName(other->getName()); int tmp; if(other->getNumberOfCells()==0) return true; @@ -1832,7 +1865,7 @@ bool MEDCouplingUMesh::areCellsIncludedIn(const MEDCouplingUMesh *other, int com * \param arr is an output parameter that returns a \b newly created instance. This array is of size 'other->getNumberOfCells()'. * \return If \a other is fully included in 'this 'true is returned. If not false is returned. */ -bool MEDCouplingUMesh::areCellsIncludedIn2(const MEDCouplingUMesh *other, DataArrayInt *& arr) const throw(INTERP_KERNEL::Exception) +bool MEDCouplingUMesh::areCellsIncludedIn2(const MEDCouplingUMesh *other, DataArrayInt *& arr) const { MEDCouplingAutoRefCountObjectPtr mesh=MergeUMeshesOnSameCoords(this,other); DataArrayInt *commonCells=0,*commonCellsI=0; @@ -1860,7 +1893,7 @@ bool MEDCouplingUMesh::areCellsIncludedIn2(const MEDCouplingUMesh *other, DataAr } } } - arr2->setName(other->getName().c_str()); + arr2->setName(other->getName()); if(arr2->presenceOfValue(0)) return false; arr=arr2.retn(); @@ -1891,7 +1924,7 @@ MEDCouplingPointSet *MEDCouplingUMesh::mergeMyselfWithOnSameCoords(const MEDCoup * \warning This method modifies can generate an unstructured mesh whose cells are not sorted by geometric type order. * In view of the MED file writing, a renumbering of cells of returned unstructured mesh (using MEDCouplingUMesh::sortCellsInMEDFileFrmt) should be necessary. */ -MEDCouplingPointSet *MEDCouplingUMesh::buildPartOfMySelf2(int start, int end, int step, bool keepCoords) const throw(INTERP_KERNEL::Exception) +MEDCouplingPointSet *MEDCouplingUMesh::buildPartOfMySelf2(int start, int end, int step, bool keepCoords) const { if(getMeshDimension()!=-1) return MEDCouplingPointSet::buildPartOfMySelf2(start,end,step,keepCoords); @@ -1955,7 +1988,7 @@ MEDCouplingPointSet *MEDCouplingUMesh::buildPartOfMySelf(const int *begin, const * \param [in] otherOnSameCoordsThanThis an another mesh with same meshdimension than \b this with exactly the same number of cells than cell ids list in [\b cellIdsBg, \b cellIdsEnd ). * Coordinate pointer of \b this and those of \b otherOnSameCoordsThanThis must be the same */ -void MEDCouplingUMesh::setPartOfMySelf(const int *cellIdsBg, const int *cellIdsEnd, const MEDCouplingUMesh& otherOnSameCoordsThanThis) throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::setPartOfMySelf(const int *cellIdsBg, const int *cellIdsEnd, const MEDCouplingUMesh& otherOnSameCoordsThanThis) { checkConnectivityFullyDefined(); otherOnSameCoordsThanThis.checkConnectivityFullyDefined(); @@ -2004,7 +2037,7 @@ void MEDCouplingUMesh::setPartOfMySelf(const int *cellIdsBg, const int *cellIdsE } } -void MEDCouplingUMesh::setPartOfMySelf2(int start, int end, int step, const MEDCouplingUMesh& otherOnSameCoordsThanThis) throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::setPartOfMySelf2(int start, int end, int step, const MEDCouplingUMesh& otherOnSameCoordsThanThis) { checkConnectivityFullyDefined(); otherOnSameCoordsThanThis.checkConnectivityFullyDefined(); @@ -2166,7 +2199,7 @@ MEDCouplingPointSet *MEDCouplingUMesh::buildBoundaryMesh(bool keepCoords) const * A cell is detected to be on boundary if it contains one or more than one face having only one father. * This method makes the assumption that \a this is fully defined (coords,connectivity). If not an exception will be thrown. */ -DataArrayInt *MEDCouplingUMesh::findCellIdsOnBoundary() const throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::findCellIdsOnBoundary() const { checkFullyDefined(); MEDCouplingAutoRefCountObjectPtr desc=DataArrayInt::New(); @@ -2218,7 +2251,7 @@ DataArrayInt *MEDCouplingUMesh::findCellIdsOnBoundary() const throw(INTERP_KERNE * \param [out] cellIdsRk1 a newly allocated array containing cells ids of s1+s2 \b into \b cellIdsRk0 subset. To get absolute ids of s1+s2 simply invoke * cellIdsRk1->transformWithIndArr(cellIdsRk0->begin(),cellIdsRk0->end()); */ -void MEDCouplingUMesh::findCellIdsLyingOn(const MEDCouplingUMesh& otherDimM1OnSameCoords, DataArrayInt *&cellIdsRk0, DataArrayInt *&cellIdsRk1) const throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::findCellIdsLyingOn(const MEDCouplingUMesh& otherDimM1OnSameCoords, DataArrayInt *&cellIdsRk0, DataArrayInt *&cellIdsRk1) const { if(getCoords()!=otherDimM1OnSameCoords.getCoords()) throw INTERP_KERNEL::Exception("MEDCouplingUMesh::findCellIdsLyingOn : coordinates pointer are not the same ! Use tryToShareSameCoords method !"); @@ -2266,7 +2299,7 @@ void MEDCouplingUMesh::findCellIdsLyingOn(const MEDCouplingUMesh& otherDimM1OnSa * * \return a newly allocated mesh lying on the same coordinates than \b this. The caller has to deal with returned mesh. */ -MEDCouplingUMesh *MEDCouplingUMesh::computeSkin() const throw(INTERP_KERNEL::Exception) +MEDCouplingUMesh *MEDCouplingUMesh::computeSkin() const { MEDCouplingAutoRefCountObjectPtr desc=DataArrayInt::New(); MEDCouplingAutoRefCountObjectPtr descIndx=DataArrayInt::New(); @@ -2297,7 +2330,7 @@ DataArrayInt *MEDCouplingUMesh::findBoundaryNodes() const return skin->computeFetchedNodeIds(); } -MEDCouplingUMesh *MEDCouplingUMesh::buildUnstructured() const throw(INTERP_KERNEL::Exception) +MEDCouplingUMesh *MEDCouplingUMesh::buildUnstructured() const { incrRef(); return const_cast(this); @@ -2374,7 +2407,7 @@ void MEDCouplingUMesh::findNodesToDuplicate(const MEDCouplingUMesh& otherDimM1On * \param [in] nodeIdsToDuplicateBg begin of node ids (included) to be duplicated in connectivity only * \param [in] nodeIdsToDuplicateEnd end of node ids (excluded) to be duplicated in connectivity only */ -void MEDCouplingUMesh::duplicateNodes(const int *nodeIdsToDuplicateBg, const int *nodeIdsToDuplicateEnd) throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::duplicateNodes(const int *nodeIdsToDuplicateBg, const int *nodeIdsToDuplicateEnd) { int nbOfNodes=getNumberOfNodes(); duplicateNodesInCoords(nodeIdsToDuplicateBg,nodeIdsToDuplicateEnd); @@ -2420,7 +2453,7 @@ void MEDCouplingUMesh::renumberNodesInConn(const int *newNodeNumbersO2N) * * \param [in] delta specifies the shift size applied to nodeId in nodal connectivity in \b this. */ -void MEDCouplingUMesh::shiftNodeNumbersInConn(int delta) throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::shiftNodeNumbersInConn(int delta) { checkConnectivityFullyDefined(); int *conn=getNodalConnectivity()->getPointer(); @@ -2455,7 +2488,7 @@ void MEDCouplingUMesh::shiftNodeNumbersInConn(int delta) throw(INTERP_KERNEL::Ex * \param [in] nodeIdsToDuplicateEnd end of node ids (excluded) to be duplicated in connectivity only * \param [in] offset the offset applied to all node ids in connectivity that are in [ \a nodeIdsToDuplicateBg, \a nodeIdsToDuplicateEnd ). */ -void MEDCouplingUMesh::duplicateNodesInConn(const int *nodeIdsToDuplicateBg, const int *nodeIdsToDuplicateEnd, int offset) throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::duplicateNodesInConn(const int *nodeIdsToDuplicateBg, const int *nodeIdsToDuplicateEnd, int offset) { checkConnectivityFullyDefined(); std::map m; @@ -2495,7 +2528,7 @@ void MEDCouplingUMesh::duplicateNodesInConn(const int *nodeIdsToDuplicateBg, con * * \param [in] old2NewBg is expected to be a dynamically allocated pointer of size at least equal to this->getNumberOfCells() */ -void MEDCouplingUMesh::renumberCells(const int *old2NewBg, bool check) throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::renumberCells(const int *old2NewBg, bool check) { checkConnectivityFullyDefined(); int nbCells=getNumberOfCells(); @@ -2675,7 +2708,7 @@ INTERP_KERNEL::NormalizedCellType MEDCouplingUMesh::getTypeOfCell(int cellId) co * \param [in] type the geometric type * \return cell ids in this having geometric type \a type. */ -DataArrayInt *MEDCouplingUMesh::giveCellsWithType(INTERP_KERNEL::NormalizedCellType type) const throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::giveCellsWithType(INTERP_KERNEL::NormalizedCellType type) const { MEDCouplingAutoRefCountObjectPtr ret=DataArrayInt::New(); @@ -2791,7 +2824,7 @@ std::string MEDCouplingUMesh::advancedRepr() const * This method returns a C++ code that is a dump of \a this. * This method will throw if this is not fully defined. */ -std::string MEDCouplingUMesh::cppRepr() const throw(INTERP_KERNEL::Exception) +std::string MEDCouplingUMesh::cppRepr() const { static const char coordsName[]="coords"; static const char connName[]="conn"; @@ -2825,12 +2858,12 @@ std::string MEDCouplingUMesh::reprConnectivityOfThis() const * This method analyzes the 3 arrays of \a this. For each the following behaviour is done : if the array is null a newly one is created * with number of tuples set to 0, if not the array is taken as this in the returned instance. */ -MEDCouplingUMesh *MEDCouplingUMesh::buildSetInstanceFromThis(int spaceDim) const throw(INTERP_KERNEL::Exception) +MEDCouplingUMesh *MEDCouplingUMesh::buildSetInstanceFromThis(int spaceDim) const { int mdim=getMeshDimension(); if(mdim<0) throw INTERP_KERNEL::Exception("MEDCouplingUMesh::buildSetInstanceFromThis : invalid mesh dimension ! Should be >= 0 !"); - MEDCouplingAutoRefCountObjectPtr ret=MEDCouplingUMesh::New(getName().c_str(),mdim); + MEDCouplingAutoRefCountObjectPtr ret=MEDCouplingUMesh::New(getName(),mdim); MEDCouplingAutoRefCountObjectPtr tmp1,tmp2; bool needToCpyCT=true; if(!_nodal_connec) @@ -2904,9 +2937,9 @@ int MEDCouplingUMesh::getNumberOfNodesInCell(int cellId) const * describing the cell types. * \throw If the coordinates array is not set. * \throw If the nodal connectivity of cells is not defined. - * \sa getAllTypes() + * \sa getAllGeoTypes() */ -std::set MEDCouplingUMesh::getTypesOfPart(const int *begin, const int *end) const throw(INTERP_KERNEL::Exception) +std::set MEDCouplingUMesh::getTypesOfPart(const int *begin, const int *end) const { checkFullyDefined(); std::set ret; @@ -2979,7 +3012,7 @@ void MEDCouplingUMesh::computeTypes() /*! * This method checks that all arrays are set. If yes nothing done if no an exception is thrown. */ -void MEDCouplingUMesh::checkFullyDefined() const throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::checkFullyDefined() const { if(!_nodal_connec_index || !_nodal_connec || !_coords) throw INTERP_KERNEL::Exception("Reverse nodal connectivity computation requires full connectivity and coordinates set in unstructured mesh."); @@ -2988,7 +3021,7 @@ void MEDCouplingUMesh::checkFullyDefined() const throw(INTERP_KERNEL::Exception) /*! * This method checks that all connectivity arrays are set. If yes nothing done if no an exception is thrown. */ -void MEDCouplingUMesh::checkConnectivityFullyDefined() const throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::checkConnectivityFullyDefined() const { if(!_nodal_connec_index || !_nodal_connec) throw INTERP_KERNEL::Exception("Reverse nodal connectivity computation requires full connectivity set in unstructured mesh."); @@ -3218,7 +3251,7 @@ MEDCouplingFieldDouble *MEDCouplingUMesh::getMeasureField(bool isAbs) const name+=getName(); int nbelem=getNumberOfCells(); MEDCouplingAutoRefCountObjectPtr field=MEDCouplingFieldDouble::New(ON_CELLS,ONE_TIME); - field->setName(name.c_str()); + field->setName(name); MEDCouplingAutoRefCountObjectPtr array=DataArrayDouble::New(); array->alloc(nbelem,1); double *area_vol=array->getPointer(); @@ -3273,7 +3306,7 @@ DataArrayDouble *MEDCouplingUMesh::getPartMeasureField(bool isAbs, const int *be name+=getName(); int nbelem=(int)std::distance(begin,end); MEDCouplingAutoRefCountObjectPtr array=DataArrayDouble::New(); - array->setName(name.c_str()); + array->setName(name); array->alloc(nbelem,1); double *area_vol=array->getPointer(); if(getMeshDimension()!=-1) @@ -3549,7 +3582,7 @@ MEDCouplingFieldDouble *MEDCouplingUMesh::buildDirectionVectorField() const * \throw If the plane does not intersect any 3D cell of \a this mesh. * \throw If \a this includes quadratic cells. */ -MEDCouplingUMesh *MEDCouplingUMesh::buildSlice3D(const double *origin, const double *vec, double eps, DataArrayInt *&cellIds) const throw(INTERP_KERNEL::Exception) +MEDCouplingUMesh *MEDCouplingUMesh::buildSlice3D(const double *origin, const double *vec, double eps, DataArrayInt *&cellIds) const { checkFullyDefined(); if(getMeshDimension()!=3 || getSpaceDimension()!=3) @@ -3615,7 +3648,7 @@ the result mesh. * \throw If the plane does not intersect any 2D cell of \a this mesh. * \throw If \a this includes quadratic cells. */ -MEDCouplingUMesh *MEDCouplingUMesh::buildSlice3DSurf(const double *origin, const double *vec, double eps, DataArrayInt *&cellIds) const throw(INTERP_KERNEL::Exception) +MEDCouplingUMesh *MEDCouplingUMesh::buildSlice3DSurf(const double *origin, const double *vec, double eps, DataArrayInt *&cellIds) const { checkFullyDefined(); if(getMeshDimension()!=2 || getSpaceDimension()!=3) @@ -3696,7 +3729,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::buildSlice3DSurf(const double *origin, const * \throw If magnitude of \a vec is less than 1e-6. * \sa buildSlice3D() */ -DataArrayInt *MEDCouplingUMesh::getCellIdsCrossingPlane(const double *origin, const double *vec, double eps) const throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::getCellIdsCrossingPlane(const double *origin, const double *vec, double eps) const { checkFullyDefined(); if(getSpaceDimension()!=3) @@ -3735,7 +3768,7 @@ DataArrayInt *MEDCouplingUMesh::getCellIdsCrossingPlane(const double *origin, co * A 1D mesh is said contiguous if : a cell i with nodal connectivity (k,p) the cell i+1 the nodal connectivity should be (p,m) * If not false is returned. In case that false is returned a call to ParaMEDMEM::MEDCouplingUMesh::mergeNodes could be usefull. */ -bool MEDCouplingUMesh::isContiguous1D() const throw(INTERP_KERNEL::Exception) +bool MEDCouplingUMesh::isContiguous1D() const { if(getMeshDimension()!=1) throw INTERP_KERNEL::Exception("MEDCouplingUMesh::isContiguous1D : this method has a sense only for 1D mesh !"); @@ -3814,7 +3847,7 @@ void MEDCouplingUMesh::project1D(const double *pt, const double *v, double eps, * dimension - 1. * \sa DataArrayDouble::distanceToTuple, MEDCouplingUMesh::distanceToPoints */ -double MEDCouplingUMesh::distanceToPoint(const double *ptBg, const double *ptEnd, int& cellId) const throw(INTERP_KERNEL::Exception) +double MEDCouplingUMesh::distanceToPoint(const double *ptBg, const double *ptEnd, int& cellId) const { int meshDim=getMeshDimension(),spaceDim=getSpaceDimension(); if(meshDim!=spaceDim-1) @@ -3853,7 +3886,7 @@ double MEDCouplingUMesh::distanceToPoint(const double *ptBg, const double *ptEnd * \throw if mesh dimension of \a this is not equal to space dimension - 1. * \sa DataArrayDouble::distanceToTuple, MEDCouplingUMesh::distanceToPoint */ -DataArrayDouble *MEDCouplingUMesh::distanceToPoints(const DataArrayDouble *pts, DataArrayInt *& cellIds) const throw(INTERP_KERNEL::Exception) +DataArrayDouble *MEDCouplingUMesh::distanceToPoints(const DataArrayDouble *pts, DataArrayInt *& cellIds) const { if(!pts) throw INTERP_KERNEL::Exception("MEDCouplingUMesh::distanceToPoints : input points pointer is NULL !"); @@ -3924,7 +3957,7 @@ DataArrayDouble *MEDCouplingUMesh::distanceToPoints(const DataArrayDouble *pts, * \param [out] cellId that corresponds to minimal distance. If the closer node is not linked to any cell in \a this -1 is returned. * \sa MEDCouplingUMesh::distanceToPoint, MEDCouplingUMesh::distanceToPoints */ -void MEDCouplingUMesh::DistanceToPoint3DSurfAlg(const double *pt, const int *cellIdsBg, const int *cellIdsEnd, const double *coords, const int *nc, const int *ncI, double& ret0, int& cellId) throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::DistanceToPoint3DSurfAlg(const double *pt, const int *cellIdsBg, const int *cellIdsEnd, const double *coords, const int *nc, const int *ncI, double& ret0, int& cellId) { cellId=-1; ret0=std::numeric_limits::max(); @@ -3963,7 +3996,7 @@ void MEDCouplingUMesh::DistanceToPoint3DSurfAlg(const double *pt, const int *cel * \param [out] cellId that corresponds to minimal distance. If the closer node is not linked to any cell in \a this -1 is returned. * \sa MEDCouplingUMesh::distanceToPoint, MEDCouplingUMesh::distanceToPoints */ -void MEDCouplingUMesh::DistanceToPoint2DCurveAlg(const double *pt, const int *cellIdsBg, const int *cellIdsEnd, const double *coords, const int *nc, const int *ncI, double& ret0, int& cellId) throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::DistanceToPoint2DCurveAlg(const double *pt, const int *cellIdsBg, const int *cellIdsEnd, const double *coords, const int *nc, const int *ncI, double& ret0, int& cellId) { cellId=-1; ret0=std::numeric_limits::max(); @@ -3988,6 +4021,9 @@ void MEDCouplingUMesh::DistanceToPoint2DCurveAlg(const double *pt, const int *ce /*! * Finds cells in contact with a ball (i.e. a point with precision). + * For speed reasons, the INTERP_KERNEL::NORM_QUAD4, INTERP_KERNEL::NORM_TRI6 and INTERP_KERNEL::NORM_QUAD8 cells are considered as convex cells to detect if a point is IN or OUT. + * If it is not the case, please change their types to INTERP_KERNEL::NORM_POLYGON or INTERP_KERNEL::NORM_QPOLYG before invoking this method. + * * \warning This method is suitable if the caller intends to evaluate only one * point, for more points getCellsContainingPoints() is recommended as it is * faster. @@ -4009,12 +4045,14 @@ int MEDCouplingUMesh::getCellContainingPoint(const double *pos, double eps) cons /*! * Finds cells in contact with a ball (i.e. a point with precision). + * For speed reasons, the INTERP_KERNEL::NORM_QUAD4, INTERP_KERNEL::NORM_TRI6 and INTERP_KERNEL::NORM_QUAD8 cells are considered as convex cells to detect if a point is IN or OUT. + * If it is not the case, please change their types to INTERP_KERNEL::NORM_POLYGON or INTERP_KERNEL::NORM_QPOLYG before invoking this method. * \warning This method is suitable if the caller intends to evaluate only one * point, for more points getCellsContainingPoints() is recommended as it is * faster. * \param [in] pos - array of coordinates of the ball central point. * \param [in] eps - ball radius. - * \param [in,out] elts - vector returning ids of the found cells. It is cleared + * \param [out] elts - vector returning ids of the found cells. It is cleared * before inserting ids. * \throw If the coordinates array is not set. * \throw If \a this->getMeshDimension() != \a this->getSpaceDimension(). @@ -4024,8 +4062,9 @@ int MEDCouplingUMesh::getCellContainingPoint(const double *pos, double eps) cons */ void MEDCouplingUMesh::getCellsContainingPoint(const double *pos, double eps, std::vector& elts) const { - std::vector eltsIndex; - getCellsContainingPoints(pos,1,eps,elts,eltsIndex); + MEDCouplingAutoRefCountObjectPtr eltsUg,eltsIndexUg; + getCellsContainingPoints(pos,1,eps,eltsUg,eltsIndexUg); + elts.clear(); elts.insert(elts.end(),eltsUg->begin(),eltsUg->end()); } /// @cond INTERNAL @@ -4064,6 +4103,7 @@ namespace ParaMEDMEM INTERP_KERNEL::EdgeLin *e1=new INTERP_KERNEL::EdgeLin(mapp2[bg[0]].first,mapp2[bg[2]].first); INTERP_KERNEL::EdgeLin *e2=new INTERP_KERNEL::EdgeLin(mapp2[bg[2]].first,mapp2[bg[1]].first); INTERP_KERNEL::SegSegIntersector inters(*e1,*e2); + // is the SEG3 degenerated, and thus can be reduced to a SEG2? bool colinearity=inters.areColinears(); delete e1; delete e2; if(colinearity) @@ -4080,11 +4120,14 @@ namespace ParaMEDMEM } /*! - * This method creates a sub mesh in Geometric2D DS. The sub mesh is composed be the sub set of cells in 'candidates' and the global mesh 'mDesc'. - * The input meth 'mDesc' must be so that mDim==1 et spaceDim==3. - * 'mapp' contains a mapping between local numbering in submesh and the global node numbering in 'mDesc'. + * This method creates a sub mesh in Geometric2D DS. The sub mesh is composed by the sub set of cells in 'candidates' taken from + * the global mesh 'mDesc'. + * The input mesh 'mDesc' must be so that mDim==1 and spaceDim==2. + * 'mapp' returns a mapping between local numbering in submesh (represented by a Node*) and the global node numbering in 'mDesc'. */ - INTERP_KERNEL::QuadraticPolygon *MEDCouplingUMeshBuildQPFromMesh(const MEDCouplingUMesh *mDesc, const std::vector& candidates, std::map& mapp) throw(INTERP_KERNEL::Exception) + INTERP_KERNEL::QuadraticPolygon *MEDCouplingUMeshBuildQPFromMesh(const MEDCouplingUMesh *mDesc, const std::vector& candidates, + std::map& mapp) + throw(INTERP_KERNEL::Exception) { mapp.clear(); std::map > mapp2;//bool is for a flag specifying if node is boundary (true) or only a middle for SEG3. @@ -4129,6 +4172,9 @@ namespace ParaMEDMEM return new INTERP_KERNEL::Node(coo1[2*nodeId],coo1[2*nodeId+1]); } + /** + * Construct a mapping between set of Nodes and the standart MEDCoupling connectivity format (c, cI). + */ void MEDCouplingUMeshBuildQPFromMesh3(const double *coo1, int offset1, const double *coo2, int offset2, const std::vector& addCoo, const int *desc1Bg, const int *desc1End, const std::vector >& intesctEdges1, /*output*/std::map& mapp, std::map& mappRev) @@ -4154,12 +4200,11 @@ namespace ParaMEDMEM template void MEDCouplingUMesh::getCellsContainingPointsAlg(const double *coords, const double *pos, int nbOfPoints, - double eps, std::vector& elts, std::vector& eltsIndex) const + double eps, MEDCouplingAutoRefCountObjectPtr& elts, MEDCouplingAutoRefCountObjectPtr& eltsIndex) const { - eltsIndex.resize(nbOfPoints+1); - eltsIndex[0]=0; - elts.clear(); - MEDCouplingAutoRefCountObjectPtr bboxArr(getBoundingBoxForBBTree()); + elts=DataArrayInt::New(); eltsIndex=DataArrayInt::New(); eltsIndex->alloc(nbOfPoints+1,1); eltsIndex->setIJ(0,0,0); elts->alloc(0,1); + int *eltsIndexPtr(eltsIndex->getPointer()); + MEDCouplingAutoRefCountObjectPtr bboxArr(getBoundingBoxForBBTree(eps)); const double *bbox(bboxArr->begin()); int nbOfCells=getNumberOfCells(); const int *conn=_nodal_connec->getConstPointer(); @@ -4168,7 +4213,7 @@ void MEDCouplingUMesh::getCellsContainingPointsAlg(const double *coords, const d BBTree myTree(&bbox[0],0,0,nbOfCells,-eps); for(int i=0;i::const_iterator iter=candidates.begin();iter!=candidates.end();iter++) { - int sz=connI[(*iter)+1]-connI[*iter]-1; - if(INTERP_KERNEL::PointLocatorAlgos >::isElementContainsPoint(pos+i*SPACEDIM, - (INTERP_KERNEL::NormalizedCellType)conn[connI[*iter]], - coords,conn+connI[*iter]+1,sz,eps)) + int sz(connI[(*iter)+1]-connI[*iter]-1); + INTERP_KERNEL::NormalizedCellType ct((INTERP_KERNEL::NormalizedCellType)conn[connI[*iter]]); + bool status(false); + if(ct!=INTERP_KERNEL::NORM_POLYGON && ct!=INTERP_KERNEL::NORM_QPOLYG) + status=INTERP_KERNEL::PointLocatorAlgos >::isElementContainsPoint(pos+i*SPACEDIM,ct,coords,conn+connI[*iter]+1,sz,eps); + else { - eltsIndex[i+1]++; - elts.push_back(*iter); + if(SPACEDIM!=2) + throw INTERP_KERNEL::Exception("MEDCouplingUMesh::getCellsContainingPointsAlg : not implemented yet for POLYGON and QPOLYGON in spaceDim 3 !"); + INTERP_KERNEL::QUADRATIC_PLANAR::_precision=eps; + INTERP_KERNEL::QUADRATIC_PLANAR::_arc_detection_precision=eps; + std::vector nodes(sz); + INTERP_KERNEL::QuadraticPolygon *pol(0); + for(int j=0;jnormalizeMe(b,c); n->applySimilarity(b,c,a); + status=pol->isInOrOut2(n); + delete pol; n->decrRef(); + } + if(status) + { + eltsIndexPtr[i+1]++; + elts->pushBackSilent(*iter); } } } @@ -4193,13 +4263,15 @@ void MEDCouplingUMesh::getCellsContainingPointsAlg(const double *coords, const d * Finds cells in contact with several balls (i.e. points with precision). * This method is an extension of getCellContainingPoint() and * getCellsContainingPoint() for the case of multiple points. + * For speed reasons, the INTERP_KERNEL::NORM_QUAD4, INTERP_KERNEL::NORM_TRI6 and INTERP_KERNEL::NORM_QUAD8 cells are considered as convex cells to detect if a point is IN or OUT. + * If it is not the case, please change their types to INTERP_KERNEL::NORM_POLYGON or INTERP_KERNEL::NORM_QPOLYG before invoking this method. * \param [in] pos - an array of coordinates of points in full interlace mode : * X0,Y0,Z0,X1,Y1,Z1,... Size of the array must be \a * this->getSpaceDimension() * \a nbOfPoints * \param [in] nbOfPoints - number of points to locate within \a this mesh. * \param [in] eps - radius of balls (i.e. the precision). - * \param [in,out] elts - vector returning ids of found cells. - * \param [in,out] eltsIndex - an array, of length \a nbOfPoints + 1, + * \param [out] elts - vector returning ids of found cells. + * \param [out] eltsIndex - an array, of length \a nbOfPoints + 1, * dividing cell ids in \a elts into groups each referring to one * point. Its every element (except the last one) is an index pointing to the * first id of a group of cells. For example cells in contact with the *i*-th @@ -4215,7 +4287,7 @@ void MEDCouplingUMesh::getCellsContainingPointsAlg(const double *coords, const d * \ref py_mcumesh_getCellsContainingPoints "Here is a Python example". */ void MEDCouplingUMesh::getCellsContainingPoints(const double *pos, int nbOfPoints, double eps, - std::vector& elts, std::vector& eltsIndex) const + MEDCouplingAutoRefCountObjectPtr& elts, MEDCouplingAutoRefCountObjectPtr& eltsIndex) const { int spaceDim=getSpaceDimension(); int mDim=getMeshDimension(); @@ -4307,7 +4379,7 @@ void MEDCouplingUMesh::checkButterflyCells(std::vector& cells, double eps) * * \return a newly allocated array containing cellIds that have been modified if any. If no cells have been impacted by this method NULL is returned. */ -DataArrayInt *MEDCouplingUMesh::convexEnvelop2D() throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::convexEnvelop2D() { if(getMeshDimension()!=2 || getSpaceDimension()!=2) throw INTERP_KERNEL::Exception("MEDCouplingUMesh::convexEnvelop2D works only for meshDim=2 and spaceDim=2 !"); @@ -4402,7 +4474,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::buildExtrudedMesh(const MEDCouplingUMesh *me * if cut3DCurve[i]==-1, it means that cell#i has been already detected to be fully part of plane defined by ('origin','vec'). * This method will throw an exception if \a this contains a non linear segment. */ -void MEDCouplingUMesh::split3DCurveWithPlane(const double *origin, const double *vec, double eps, std::vector& cut3DCurve) throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::split3DCurveWithPlane(const double *origin, const double *vec, double eps, std::vector& cut3DCurve) { checkFullyDefined(); if(getMeshDimension()!=1 || getSpaceDimension()!=3) @@ -4509,7 +4581,7 @@ DataArrayDouble *MEDCouplingUMesh::fillExtCoordsUsingTranslation(const MEDCoupli * \param mesh1D is the input 1D mesh used for translation and automatic rotation computation. * \return newCoords new coords filled by this method. */ -DataArrayDouble *MEDCouplingUMesh::fillExtCoordsUsingTranslAndAutoRotation(const MEDCouplingUMesh *mesh1D, bool isQuad) const throw(INTERP_KERNEL::Exception) +DataArrayDouble *MEDCouplingUMesh::fillExtCoordsUsingTranslAndAutoRotation(const MEDCouplingUMesh *mesh1D, bool isQuad) const { if(mesh1D->getSpaceDimension()==2) return fillExtCoordsUsingTranslAndAutoRotation2D(mesh1D,isQuad); @@ -4523,7 +4595,7 @@ DataArrayDouble *MEDCouplingUMesh::fillExtCoordsUsingTranslAndAutoRotation(const * \param mesh1D is the input 1D mesh used for translation and automatic rotation computation. * \return newCoords new coords filled by this method. */ -DataArrayDouble *MEDCouplingUMesh::fillExtCoordsUsingTranslAndAutoRotation2D(const MEDCouplingUMesh *mesh1D, bool isQuad) const throw(INTERP_KERNEL::Exception) +DataArrayDouble *MEDCouplingUMesh::fillExtCoordsUsingTranslAndAutoRotation2D(const MEDCouplingUMesh *mesh1D, bool isQuad) const { if(isQuad) throw INTERP_KERNEL::Exception("MEDCouplingUMesh::fillExtCoordsUsingTranslAndAutoRotation2D : not implemented for quadratic cells !"); @@ -4567,7 +4639,7 @@ DataArrayDouble *MEDCouplingUMesh::fillExtCoordsUsingTranslAndAutoRotation2D(con * \param mesh1D is the input 1D mesh used for translation and automatic rotation computation. * \return newCoords new coords filled by this method. */ -DataArrayDouble *MEDCouplingUMesh::fillExtCoordsUsingTranslAndAutoRotation3D(const MEDCouplingUMesh *mesh1D, bool isQuad) const throw(INTERP_KERNEL::Exception) +DataArrayDouble *MEDCouplingUMesh::fillExtCoordsUsingTranslAndAutoRotation3D(const MEDCouplingUMesh *mesh1D, bool isQuad) const { if(isQuad) throw INTERP_KERNEL::Exception("MEDCouplingUMesh::fillExtCoordsUsingTranslAndAutoRotation3D : not implemented for quadratic cells !"); @@ -4727,7 +4799,7 @@ bool MEDCouplingUMesh::isPresenceOfQuadratic() const * \throw If the coordinates array is not set. * \throw If the nodal connectivity of cells is not defined. */ -void MEDCouplingUMesh::convertQuadraticCellsToLinear() throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::convertQuadraticCellsToLinear() { checkFullyDefined(); int nbOfCells=getNumberOfCells(); @@ -4800,7 +4872,7 @@ void MEDCouplingUMesh::convertQuadraticCellsToLinear() throw(INTERP_KERNEL::Exce * * \sa MEDCouplingUMesh::convertQuadraticCellsToLinear */ -DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic(int conversionType) throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic(int conversionType) { DataArrayInt *conn=0,*connI=0; DataArrayDouble *coords=0; @@ -4865,7 +4937,7 @@ DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic(int conversionType * \return a newly created DataArrayInt instance that the caller should deal with containing cell ids of converted cells. * \sa MEDCouplingUMesh::convertLinearCellsToQuadratic. */ -DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic1D0(DataArrayInt *&conn, DataArrayInt *&connI, DataArrayDouble *& coords, std::set& types) const throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic1D0(DataArrayInt *&conn, DataArrayInt *&connI, DataArrayDouble *& coords, std::set& types) const { MEDCouplingAutoRefCountObjectPtr bary=getBarycenterAndOwner(); MEDCouplingAutoRefCountObjectPtr newConn=DataArrayInt::New(); newConn->alloc(0,1); @@ -4902,7 +4974,7 @@ DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic1D0(DataArrayInt *& return ret.retn(); } -DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic2DAnd3D0(const MEDCouplingUMesh *m1D, const DataArrayInt *desc, const DataArrayInt *descI, DataArrayInt *&conn, DataArrayInt *&connI, DataArrayDouble *& coords, std::set& types) const throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic2DAnd3D0(const MEDCouplingUMesh *m1D, const DataArrayInt *desc, const DataArrayInt *descI, DataArrayInt *&conn, DataArrayInt *&connI, DataArrayDouble *& coords, std::set& types) const { MEDCouplingAutoRefCountObjectPtr newConn=DataArrayInt::New(); newConn->alloc(0,1); MEDCouplingAutoRefCountObjectPtr newConnI=DataArrayInt::New(); newConnI->alloc(1,1); newConnI->setIJ(0,0,0); @@ -4953,7 +5025,7 @@ DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic2DAnd3D0(const MEDC * \return a newly created DataArrayInt instance that the caller should deal with containing cell ids of converted cells. * \sa MEDCouplingUMesh::convertLinearCellsToQuadratic. */ -DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic2D0(DataArrayInt *&conn, DataArrayInt *&connI, DataArrayDouble *& coords, std::set& types) const throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic2D0(DataArrayInt *&conn, DataArrayInt *&connI, DataArrayDouble *& coords, std::set& types) const { MEDCouplingAutoRefCountObjectPtr desc(DataArrayInt::New()),descI(DataArrayInt::New()),tmp2(DataArrayInt::New()),tmp3(DataArrayInt::New()); @@ -4961,7 +5033,7 @@ DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic2D0(DataArrayInt *& return convertLinearCellsToQuadratic2DAnd3D0(m1D,desc,descI,conn,connI,coords,types); } -DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic2D1(DataArrayInt *&conn, DataArrayInt *&connI, DataArrayDouble *& coords, std::set& types) const throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic2D1(DataArrayInt *&conn, DataArrayInt *&connI, DataArrayDouble *& coords, std::set& types) const { MEDCouplingAutoRefCountObjectPtr desc(DataArrayInt::New()),descI(DataArrayInt::New()),tmp2(DataArrayInt::New()),tmp3(DataArrayInt::New()); MEDCouplingAutoRefCountObjectPtr m1D=buildDescendingConnectivity(desc,descI,tmp2,tmp3); tmp2=0; tmp3=0; @@ -5018,14 +5090,14 @@ DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic2D1(DataArrayInt *& * \return a newly created DataArrayInt instance that the caller should deal with containing cell ids of converted cells. * \sa MEDCouplingUMesh::convertLinearCellsToQuadratic. */ -DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic3D0(DataArrayInt *&conn, DataArrayInt *&connI, DataArrayDouble *& coords, std::set& types) const throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic3D0(DataArrayInt *&conn, DataArrayInt *&connI, DataArrayDouble *& coords, std::set& types) const { MEDCouplingAutoRefCountObjectPtr desc(DataArrayInt::New()),descI(DataArrayInt::New()),tmp2(DataArrayInt::New()),tmp3(DataArrayInt::New()); MEDCouplingAutoRefCountObjectPtr m1D=explode3DMeshTo1D(desc,descI,tmp2,tmp3); tmp2=0; tmp3=0; return convertLinearCellsToQuadratic2DAnd3D0(m1D,desc,descI,conn,connI,coords,types); } -DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic3D1(DataArrayInt *&conn, DataArrayInt *&connI, DataArrayDouble *& coords, std::set& types) const throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic3D1(DataArrayInt *&conn, DataArrayInt *&connI, DataArrayDouble *& coords, std::set& types) const { MEDCouplingAutoRefCountObjectPtr desc2(DataArrayInt::New()),desc2I(DataArrayInt::New()),tmp2(DataArrayInt::New()),tmp3(DataArrayInt::New()); MEDCouplingAutoRefCountObjectPtr m2D=buildDescendingConnectivityGen(desc2,desc2I,tmp2,tmp3,MEDCouplingFastNbrer); tmp2=0; tmp3=0; @@ -5118,7 +5190,7 @@ DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic3D1(DataArrayInt *& * \throw If \a this->getMeshDimension() != 2. * \throw If \a this->getSpaceDimension() != 2. */ -void MEDCouplingUMesh::tessellate2D(double eps) throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::tessellate2D(double eps) { checkFullyDefined(); if(getMeshDimension()!=2 || getSpaceDimension()!=2) @@ -5147,7 +5219,7 @@ void MEDCouplingUMesh::tessellate2D(double eps) throw(INTERP_KERNEL::Exception) * \throw If \a this->getMeshDimension() != 1. * \throw If \a this->getSpaceDimension() != 2. */ -void MEDCouplingUMesh::tessellate2DCurve(double eps) throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::tessellate2DCurve(double eps) { checkFullyDefined(); if(getMeshDimension()!=1 || getSpaceDimension()!=2) @@ -5202,7 +5274,7 @@ void MEDCouplingUMesh::tessellate2DCurve(double eps) throw(INTERP_KERNEL::Except newConnIPtr[1]=newConnIPtr[0]+3; } } - if(addCoo.empty() && ((int)newConn.size())==_nodal_connec->getNumberOfTuples())//nothing happens during tasselation : no update needed + if(addCoo.empty() && ((int)newConn.size())==_nodal_connec->getNumberOfTuples())//nothing happens during tessellation : no update needed return ; _types=types; DataArrayInt::SetArrayIn(newConnI,_nodal_connec_index); @@ -5239,9 +5311,9 @@ 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 + * \sa MEDCouplingUMesh::tetrahedrize, MEDCoupling1SGTUMesh::sortHexa8EachOther */ -DataArrayInt *MEDCouplingUMesh::simplexize(int policy) throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::simplexize(int policy) { switch(policy) { @@ -5272,7 +5344,7 @@ DataArrayInt *MEDCouplingUMesh::simplexize(int policy) throw(INTERP_KERNEL::Exce * \throw If the nodal connectivity of cells is not defined. * \throw If \a this->getMeshDimension() < 1. */ -bool MEDCouplingUMesh::areOnlySimplexCells() const throw(INTERP_KERNEL::Exception) +bool MEDCouplingUMesh::areOnlySimplexCells() const { checkFullyDefined(); int mdim=getMeshDimension(); @@ -5293,7 +5365,7 @@ bool MEDCouplingUMesh::areOnlySimplexCells() const throw(INTERP_KERNEL::Exceptio /*! * This method implements policy 0 of virtual method ParaMEDMEM::MEDCouplingUMesh::simplexize. */ -DataArrayInt *MEDCouplingUMesh::simplexizePol0() throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::simplexizePol0() { checkConnectivityFullyDefined(); if(getMeshDimension()!=2) @@ -5346,7 +5418,7 @@ DataArrayInt *MEDCouplingUMesh::simplexizePol0() throw(INTERP_KERNEL::Exception) /*! * This method implements policy 1 of virtual method ParaMEDMEM::MEDCouplingUMesh::simplexize. */ -DataArrayInt *MEDCouplingUMesh::simplexizePol1() throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::simplexizePol1() { checkConnectivityFullyDefined(); if(getMeshDimension()!=2) @@ -5399,7 +5471,7 @@ DataArrayInt *MEDCouplingUMesh::simplexizePol1() throw(INTERP_KERNEL::Exception) /*! * This method implements policy INTERP_KERNEL::PLANAR_FACE_5 of virtual method ParaMEDMEM::MEDCouplingUMesh::simplexize. */ -DataArrayInt *MEDCouplingUMesh::simplexizePlanarFace5() throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::simplexizePlanarFace5() { checkConnectivityFullyDefined(); if(getMeshDimension()!=3) @@ -5451,7 +5523,7 @@ DataArrayInt *MEDCouplingUMesh::simplexizePlanarFace5() throw(INTERP_KERNEL::Exc /*! * This method implements policy INTERP_KERNEL::PLANAR_FACE_6 of virtual method ParaMEDMEM::MEDCouplingUMesh::simplexize. */ -DataArrayInt *MEDCouplingUMesh::simplexizePlanarFace6() throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::simplexizePlanarFace6() { checkConnectivityFullyDefined(); if(getMeshDimension()!=3) @@ -5508,7 +5580,7 @@ DataArrayInt *MEDCouplingUMesh::simplexizePlanarFace6() throw(INTERP_KERNEL::Exc * \param desc is descending connectivity in format specified in MEDCouplingUMesh::buildDescendingConnectivity2 * \param descIndex is descending connectivity index in format specified in MEDCouplingUMesh::buildDescendingConnectivity2 */ -void MEDCouplingUMesh::subDivide2DMesh(const int *nodeSubdived, const int *nodeIndxSubdived, const int *desc, const int *descIndex) throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::subDivide2DMesh(const int *nodeSubdived, const int *nodeIndxSubdived, const int *desc, const int *descIndex) { checkFullyDefined(); if(getMeshDimension()!=2) @@ -5592,7 +5664,7 @@ void MEDCouplingUMesh::subDivide2DMesh(const int *nodeSubdived, const int *nodeI * \throw If the coordinates array is not set. * \throw If the nodal connectivity of cells is not defined. */ -void MEDCouplingUMesh::convertDegeneratedCells() throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::convertDegeneratedCells() { checkFullyDefined(); if(getMeshDimension()<=1) @@ -5639,7 +5711,7 @@ void MEDCouplingUMesh::convertDegeneratedCells() throw(INTERP_KERNEL::Exception) * \ref cpp_mcumesh_are2DCellsNotCorrectlyOriented "Here is a C++ example".
* \ref py_mcumesh_are2DCellsNotCorrectlyOriented "Here is a Python example". */ -void MEDCouplingUMesh::are2DCellsNotCorrectlyOriented(const double *vec, bool polyOnly, std::vector& cells) const throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::are2DCellsNotCorrectlyOriented(const double *vec, bool polyOnly, std::vector& cells) const { if(getMeshDimension()!=2 || getSpaceDimension()!=3) throw INTERP_KERNEL::Exception("Invalid mesh to apply are2DCellsNotCorrectlyOriented on it : must be meshDim==2 and spaceDim==3 !"); @@ -5673,7 +5745,7 @@ void MEDCouplingUMesh::are2DCellsNotCorrectlyOriented(const double *vec, bool po * \ref cpp_mcumesh_are2DCellsNotCorrectlyOriented "Here is a C++ example".
* \ref py_mcumesh_are2DCellsNotCorrectlyOriented "Here is a Python example". */ -void MEDCouplingUMesh::orientCorrectly2DCells(const double *vec, bool polyOnly) throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::orientCorrectly2DCells(const double *vec, bool polyOnly) { if(getMeshDimension()!=2 || getSpaceDimension()!=3) throw INTERP_KERNEL::Exception("Invalid mesh to apply orientCorrectly2DCells on it : must be meshDim==2 and spaceDim==3 !"); @@ -5687,13 +5759,25 @@ void MEDCouplingUMesh::orientCorrectly2DCells(const double *vec, bool polyOnly) INTERP_KERNEL::NormalizedCellType type=(INTERP_KERNEL::NormalizedCellType)conn[connI[i]]; if(!polyOnly || (type==INTERP_KERNEL::NORM_POLYGON || type==INTERP_KERNEL::NORM_QPOLYG)) { - bool isQuadratic=INTERP_KERNEL::CellModel::GetCellModel(type).isQuadratic(); + bool isQuadratic(INTERP_KERNEL::CellModel::GetCellModel(type).isQuadratic()); if(!IsPolygonWellOriented(isQuadratic,vec,conn+connI[i]+1,conn+connI[i+1],coordsPtr)) { isModified=true; - std::vector tmp(connI[i+1]-connI[i]-2); - std::copy(conn+connI[i]+2,conn+connI[i+1],tmp.rbegin()); - std::copy(tmp.begin(),tmp.end(),conn+connI[i]+2); + if(!isQuadratic) + { + std::vector tmp(connI[i+1]-connI[i]-2); + std::copy(conn+connI[i]+2,conn+connI[i+1],tmp.rbegin()); + std::copy(tmp.begin(),tmp.end(),conn+connI[i]+2); + } + else + { + int sz(((int)(connI[i+1]-connI[i]-1))/2); + std::vector tmp0(sz-1),tmp1(sz); + std::copy(conn+connI[i]+2,conn+connI[i]+1+sz,tmp0.rbegin()); + std::copy(conn+connI[i]+1+sz,conn+connI[i+1],tmp1.rbegin()); + std::copy(tmp0.begin(),tmp0.end(),conn+connI[i]+2); + std::copy(tmp1.begin(),tmp1.end(),conn+connI[i]+1+sz); + } } } } @@ -5715,7 +5799,7 @@ void MEDCouplingUMesh::orientCorrectly2DCells(const double *vec, bool polyOnly) * \ref cpp_mcumesh_arePolyhedronsNotCorrectlyOriented "Here is a C++ example".
* \ref py_mcumesh_arePolyhedronsNotCorrectlyOriented "Here is a Python example". */ -void MEDCouplingUMesh::arePolyhedronsNotCorrectlyOriented(std::vector& cells) const throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::arePolyhedronsNotCorrectlyOriented(std::vector& cells) const { if(getMeshDimension()!=3 || getSpaceDimension()!=3) throw INTERP_KERNEL::Exception("Invalid mesh to apply arePolyhedronsNotCorrectlyOriented on it : must be meshDim==3 and spaceDim==3 !"); @@ -5747,7 +5831,7 @@ void MEDCouplingUMesh::arePolyhedronsNotCorrectlyOriented(std::vector& cell * \ref py_mcumesh_arePolyhedronsNotCorrectlyOriented "Here is a Python example". * \sa MEDCouplingUMesh::findAndCorrectBadOriented3DCells */ -void MEDCouplingUMesh::orientCorrectlyPolyhedrons() throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::orientCorrectlyPolyhedrons() { if(getMeshDimension()!=3 || getSpaceDimension()!=3) throw INTERP_KERNEL::Exception("Invalid mesh to apply orientCorrectlyPolyhedrons on it : must be meshDim==3 and spaceDim==3 !"); @@ -5792,7 +5876,7 @@ void MEDCouplingUMesh::orientCorrectlyPolyhedrons() throw(INTERP_KERNEL::Excepti * \ref py_mcumesh_findAndCorrectBadOriented3DExtrudedCells "Here is a Python example". * \sa MEDCouplingUMesh::findAndCorrectBadOriented3DCells */ -DataArrayInt *MEDCouplingUMesh::findAndCorrectBadOriented3DExtrudedCells() throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::findAndCorrectBadOriented3DExtrudedCells() { const char msg[]="check3DCellsWellOriented detection works only for 3D cells !"; if(getMeshDimension()!=3) @@ -5829,7 +5913,7 @@ DataArrayInt *MEDCouplingUMesh::findAndCorrectBadOriented3DExtrudedCells() throw * \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) +DataArrayInt *MEDCouplingUMesh::findAndCorrectBadOriented3DCells() { if(getMeshDimension()!=3 || getSpaceDimension()!=3) throw INTERP_KERNEL::Exception("Invalid mesh to apply findAndCorrectBadOriented3DCells on it : must be meshDim==3 and spaceDim==3 !"); @@ -5896,7 +5980,7 @@ DataArrayInt *MEDCouplingUMesh::findAndCorrectBadOriented3DCells() throw(INTERP_ * \param vec output of size at least 3 used to store the normal vector (with norm equal to Area ) of searched plane. * \param pos output of size at least 3 used to store a point owned of searched plane. */ -void MEDCouplingUMesh::getFastAveragePlaneOfThis(double *vec, double *pos) const throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::getFastAveragePlaneOfThis(double *vec, double *pos) const { if(getMeshDimension()!=2 || getSpaceDimension()!=3) throw INTERP_KERNEL::Exception("Invalid mesh to apply getFastAveragePlaneOfThis on it : must be meshDim==2 and spaceDim==3 !"); @@ -5930,7 +6014,7 @@ void MEDCouplingUMesh::getFastAveragePlaneOfThis(double *vec, double *pos) const * \throw If \a this->getSpaceDimension() is neither 2 nor 3. * \throw If \a this mesh includes cells of type different from the ones enumerated above. */ -MEDCouplingFieldDouble *MEDCouplingUMesh::getEdgeRatioField() const throw(INTERP_KERNEL::Exception) +MEDCouplingFieldDouble *MEDCouplingUMesh::getEdgeRatioField() const { checkCoherency(); int spaceDim=getSpaceDimension(); @@ -6002,7 +6086,7 @@ MEDCouplingFieldDouble *MEDCouplingUMesh::getEdgeRatioField() const throw(INTERP * \throw If \a this->getSpaceDimension() is neither 2 nor 3. * \throw If \a this mesh includes cells of type different from the ones enumerated above. */ -MEDCouplingFieldDouble *MEDCouplingUMesh::getAspectRatioField() const throw(INTERP_KERNEL::Exception) +MEDCouplingFieldDouble *MEDCouplingUMesh::getAspectRatioField() const { checkCoherency(); int spaceDim=getSpaceDimension(); @@ -6073,7 +6157,7 @@ MEDCouplingFieldDouble *MEDCouplingUMesh::getAspectRatioField() const throw(INTE * \throw If \a this->getSpaceDimension() != 3. * \throw If \a this mesh includes cells of type different from the ones enumerated above. */ -MEDCouplingFieldDouble *MEDCouplingUMesh::getWarpField() const throw(INTERP_KERNEL::Exception) +MEDCouplingFieldDouble *MEDCouplingUMesh::getWarpField() const { checkCoherency(); int spaceDim=getSpaceDimension(); @@ -6133,7 +6217,7 @@ MEDCouplingFieldDouble *MEDCouplingUMesh::getWarpField() const throw(INTERP_KERN * \throw If \a this->getSpaceDimension() != 3. * \throw If \a this mesh includes cells of type different from the ones enumerated above. */ -MEDCouplingFieldDouble *MEDCouplingUMesh::getSkewField() const throw(INTERP_KERNEL::Exception) +MEDCouplingFieldDouble *MEDCouplingUMesh::getSkewField() const { checkCoherency(); int spaceDim=getSpaceDimension(); @@ -6177,12 +6261,45 @@ MEDCouplingFieldDouble *MEDCouplingUMesh::getSkewField() const throw(INTERP_KERN /*! * This method aggregate the bbox of each cell and put it into bbox parameter. * - * \return DataArrayDouble * - having \a this number of cells tuples and 2*spacedim components. + * \param [in] arcDetEps - a parameter specifying in case of 2D quadratic polygon cell the detection limit between linear and arc circle. (By default 1e-12) + * For all other cases this input parameter is ignored. + * \return DataArrayDouble * - newly created object (to be managed by the caller) \a this number of cells tuples and 2*spacedim components. * * \throw If \a this is not fully set (coordinates and connectivity). * \throw If a cell in \a this has no valid nodeId. + * \sa MEDCouplingUMesh::getBoundingBoxForBBTreeFast, MEDCouplingUMesh::getBoundingBoxForBBTree2DQuadratic */ -DataArrayDouble *MEDCouplingUMesh::getBoundingBoxForBBTree() const +DataArrayDouble *MEDCouplingUMesh::getBoundingBoxForBBTree(double arcDetEps) const +{ + int mDim(getMeshDimension()),sDim(getSpaceDimension()); + if(mDim!=2 || sDim!=2) + return getBoundingBoxForBBTreeFast(); + else + { + bool presenceOfQuadratic(false); + for(std::set::const_iterator it=_types.begin();it!=_types.end();it++) + { + const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel(*it)); + if(cm.isQuadratic()) + presenceOfQuadratic=true; + } + if(presenceOfQuadratic) + return getBoundingBoxForBBTree2DQuadratic(arcDetEps); + else + return getBoundingBoxForBBTreeFast(); + } +} + +/*! + * This method aggregate the bbox of each cell only considering the nodes constituting each cell and put it into bbox parameter. + * So meshes having quadratic cells the computed bounding boxes can be invalid ! + * + * \return DataArrayDouble * - newly created object (to be managed by the caller) \a this number of cells tuples and 2*spacedim components. + * + * \throw If \a this is not fully set (coordinates and connectivity). + * \throw If a cell in \a this has no valid nodeId. + */ +DataArrayDouble *MEDCouplingUMesh::getBoundingBoxForBBTreeFast() const { checkFullyDefined(); int spaceDim(getSpaceDimension()),nbOfCells(getNumberOfCells()),nbOfNodes(getNumberOfNodes()); @@ -6198,11 +6315,10 @@ DataArrayDouble *MEDCouplingUMesh::getBoundingBoxForBBTree() const for(int i=0;i=0 && nodeId ret(DataArrayDouble::New()); ret->alloc(nbOfCells,2*spaceDim); + double *bbox(ret->getPointer()); + const double *coords(_coords->getConstPointer()); + const int *conn(_nodal_connec->getConstPointer()),*connI(_nodal_connec_index->getConstPointer()); + for(int i=0;i nodes(sz); + INTERP_KERNEL::QuadraticPolygon *pol(0); + for(int j=0;jfillBounds(b); delete pol; + bbox[0]=b.getXMin(); bbox[1]=b.getXMax(); bbox[2]=b.getYMin(); bbox[3]=b.getYMax(); + } + return ret.retn(); +} + /// @cond INTERNAL namespace ParaMEDMEMImpl @@ -6257,14 +6415,14 @@ namespace ParaMEDMEMImpl * For every k in [0,n] ret[3*k+2]==-1 because it has no sense here. * This parameter is kept only for compatibility with other methode listed above. */ -std::vector MEDCouplingUMesh::getDistributionOfTypes() const throw(INTERP_KERNEL::Exception) +std::vector MEDCouplingUMesh::getDistributionOfTypes() const { checkConnectivityFullyDefined(); const int *conn=_nodal_connec->getConstPointer(); const int *connI=_nodal_connec_index->getConstPointer(); const int *work=connI; int nbOfCells=getNumberOfCells(); - std::size_t n=getAllTypes().size(); + std::size_t n=getAllGeoTypes().size(); std::vector ret(3*n,-1); //ret[3*k+2]==-1 because it has no sense here std::set types; for(std::size_t i=0;work!=connI+nbOfCells;i++) @@ -6302,7 +6460,7 @@ std::vector MEDCouplingUMesh::getDistributionOfTypes() const throw(INTERP_K * If it exists a geometric type in \a this \b not in \a code \b no exception is thrown * and a DataArrayInt instance is returned that the user has the responsability to deallocate. */ -DataArrayInt *MEDCouplingUMesh::checkTypeConsistencyAndContig(const std::vector& code, const std::vector& idsPerType) const throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::checkTypeConsistencyAndContig(const std::vector& code, const std::vector& idsPerType) const { if(code.empty()) throw INTERP_KERNEL::Exception("MEDCouplingUMesh::checkTypeConsistencyAndContig : code is empty, should not !"); @@ -6402,7 +6560,7 @@ DataArrayInt *MEDCouplingUMesh::checkTypeConsistencyAndContig(const std::vector< * This vector can be empty in case of all geometric type cells are fully covered in ascending in the given input \a profile. * \throw if \a profile has not exactly one component. It throws too, if \a profile contains some values not in [0,getNumberOfCells()) or if \a this is not fully defined */ -void MEDCouplingUMesh::splitProfilePerType(const DataArrayInt *profile, std::vector& code, std::vector& idsInPflPerType, std::vector& idsPerType) const throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::splitProfilePerType(const DataArrayInt *profile, std::vector& code, std::vector& idsInPflPerType, std::vector& idsPerType) const { if(!profile) throw INTERP_KERNEL::Exception("MEDCouplingUMesh::splitProfilePerType : input profile is NULL !"); @@ -6479,7 +6637,7 @@ void MEDCouplingUMesh::splitProfilePerType(const DataArrayInt *profile, std::vec * The following equality should be verified 'nM1LevMesh->getMeshDimension()==this->getMeshDimension()-1' * This method returns 5+2 elements. 'desc', 'descIndx', 'revDesc', 'revDescIndx' and 'meshnM1' behaves exactly as ParaMEDMEM::MEDCouplingUMesh::buildDescendingConnectivity except the content as described after. The returned array specifies the n-1 mesh reordered by type as MEDMEM does. 'nM1LevMeshIds' contains the ids in returned 'meshnM1'. Finally 'meshnM1Old2New' contains numbering old2new that is to say the cell #k in coarse 'nM1LevMesh' will have the number ret[k] in returned mesh 'nM1LevMesh' MEDMEM reordered. */ -MEDCouplingUMesh *MEDCouplingUMesh::emulateMEDMEMBDC(const MEDCouplingUMesh *nM1LevMesh, DataArrayInt *desc, DataArrayInt *descIndx, DataArrayInt *&revDesc, DataArrayInt *&revDescIndx, DataArrayInt *& nM1LevMeshIds, DataArrayInt *&meshnM1Old2New) const throw(INTERP_KERNEL::Exception) +MEDCouplingUMesh *MEDCouplingUMesh::emulateMEDMEMBDC(const MEDCouplingUMesh *nM1LevMesh, DataArrayInt *desc, DataArrayInt *descIndx, DataArrayInt *&revDesc, DataArrayInt *&revDescIndx, DataArrayInt *& nM1LevMeshIds, DataArrayInt *&meshnM1Old2New) const { checkFullyDefined(); nM1LevMesh->checkFullyDefined(); @@ -6524,7 +6682,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::emulateMEDMEMBDC(const MEDCouplingUMesh *nM1 * this array using decrRef() as it is no more needed. * \throw If the nodal connectivity of cells is not defined. */ -DataArrayInt *MEDCouplingUMesh::sortCellsInMEDFileFrmt() throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::sortCellsInMEDFileFrmt() { checkConnectivityFullyDefined(); MEDCouplingAutoRefCountObjectPtr ret=getRenumArrForMEDFileFrmt(); @@ -6560,7 +6718,7 @@ bool MEDCouplingUMesh::checkConsecutiveCellTypes() const * * \sa MEDCouplingUMesh::checkConsecutiveCellTypesAndOrder */ -bool MEDCouplingUMesh::checkConsecutiveCellTypesForMEDFileFrmt() const throw(INTERP_KERNEL::Exception) +bool MEDCouplingUMesh::checkConsecutiveCellTypesForMEDFileFrmt() const { return checkConsecutiveCellTypesAndOrder(MEDMEM_ORDER,MEDMEM_ORDER+N_MEDMEM_ORDER); } @@ -6612,7 +6770,7 @@ bool MEDCouplingUMesh::checkConsecutiveCellTypesAndOrder(const INTERP_KERNEL::No * that tells for each cell the pos of its type in the array on type given in input parameter. The 2nd output parameter is an array with the same * number of tuples than input type array and with one component. This 2nd output array gives type by type the number of occurence of type in 'this'. */ -DataArrayInt *MEDCouplingUMesh::getLevArrPerCellTypes(const INTERP_KERNEL::NormalizedCellType *orderBg, const INTERP_KERNEL::NormalizedCellType *orderEnd, DataArrayInt *&nbPerType) const throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::getLevArrPerCellTypes(const INTERP_KERNEL::NormalizedCellType *orderBg, const INTERP_KERNEL::NormalizedCellType *orderEnd, DataArrayInt *&nbPerType) const { checkConnectivityFullyDefined(); int nbOfCells=getNumberOfCells(); @@ -6653,7 +6811,7 @@ DataArrayInt *MEDCouplingUMesh::getLevArrPerCellTypes(const INTERP_KERNEL::Norma * * \sa MEDCouplingUMesh::getRenumArrForConsecutiveCellTypesSpec, MEDCouplingUMesh::sortCellsInMEDFileFrmt. */ -DataArrayInt *MEDCouplingUMesh::getRenumArrForMEDFileFrmt() const throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::getRenumArrForMEDFileFrmt() const { return getRenumArrForConsecutiveCellTypesSpec(MEDMEM_ORDER,MEDMEM_ORDER+N_MEDMEM_ORDER); } @@ -6664,7 +6822,7 @@ DataArrayInt *MEDCouplingUMesh::getRenumArrForMEDFileFrmt() const throw(INTERP_K * The mesh after this call to MEDCouplingMesh::renumberCells will pass the test of MEDCouplingUMesh::checkConsecutiveCellTypesAndOrder with the same inputs. * The returned array minimizes the permutations that is to say the order of cells inside same geometric type remains the same. */ -DataArrayInt *MEDCouplingUMesh::getRenumArrForConsecutiveCellTypesSpec(const INTERP_KERNEL::NormalizedCellType *orderBg, const INTERP_KERNEL::NormalizedCellType *orderEnd) const throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::getRenumArrForConsecutiveCellTypesSpec(const INTERP_KERNEL::NormalizedCellType *orderBg, const INTERP_KERNEL::NormalizedCellType *orderEnd) const { DataArrayInt *nbPerType=0; MEDCouplingAutoRefCountObjectPtr tmpa=getLevArrPerCellTypes(orderBg,orderEnd,nbPerType); @@ -6749,13 +6907,13 @@ std::vector MEDCouplingUMesh::splitByType() const * \throw If the nodal connectivity of \a this is not fully defined. * \throw If the internal data is not coherent. */ -MEDCoupling1GTUMesh *MEDCouplingUMesh::convertIntoSingleGeoTypeMesh() const throw(INTERP_KERNEL::Exception) +MEDCoupling1GTUMesh *MEDCouplingUMesh::convertIntoSingleGeoTypeMesh() const { checkConnectivityFullyDefined(); if(_types.size()!=1) throw INTERP_KERNEL::Exception("MEDCouplingUMesh::convertIntoSingleGeoTypeMesh : current mesh does not contain exactly one geometric type !"); INTERP_KERNEL::NormalizedCellType typ=*_types.begin(); - MEDCouplingAutoRefCountObjectPtr ret=MEDCoupling1GTUMesh::New(getName().c_str(),typ); + MEDCouplingAutoRefCountObjectPtr ret=MEDCoupling1GTUMesh::New(getName(),typ); ret->setCoords(getCoords()); MEDCoupling1SGTUMesh *retC=dynamic_cast((MEDCoupling1GTUMesh*)ret); if(retC) @@ -6776,7 +6934,7 @@ MEDCoupling1GTUMesh *MEDCouplingUMesh::convertIntoSingleGeoTypeMesh() const thro return ret.retn(); } -DataArrayInt *MEDCouplingUMesh::convertNodalConnectivityToStaticGeoTypeMesh() const throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::convertNodalConnectivityToStaticGeoTypeMesh() const { checkConnectivityFullyDefined(); if(_types.size()!=1) @@ -6810,7 +6968,7 @@ DataArrayInt *MEDCouplingUMesh::convertNodalConnectivityToStaticGeoTypeMesh() co return connOut.retn(); } -void MEDCouplingUMesh::convertNodalConnectivityToDynamicGeoTypeMesh(DataArrayInt *&nodalConn, DataArrayInt *&nodalConnIndex) const throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::convertNodalConnectivityToDynamicGeoTypeMesh(DataArrayInt *&nodalConn, DataArrayInt *&nodalConnIndex) const { static const char msg0[]="MEDCouplingUMesh::convertNodalConnectivityToDynamicGeoTypeMesh : nodal connectivity in this are invalid ! Call checkCoherency2 !"; checkConnectivityFullyDefined(); @@ -6913,7 +7071,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::AggregateSortedByTypeMeshesOnSameCoords(cons * This method returns a newly created DataArrayInt instance. * This method retrieves cell ids in [ \a begin, \a end ) that have the type \a type. */ -DataArrayInt *MEDCouplingUMesh::keepCellIdsByType(INTERP_KERNEL::NormalizedCellType type, const int *begin, const int *end) const throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::keepCellIdsByType(INTERP_KERNEL::NormalizedCellType type, const int *begin, const int *end) const { checkFullyDefined(); const int *conn=_nodal_connec->getConstPointer(); @@ -6929,13 +7087,13 @@ DataArrayInt *MEDCouplingUMesh::keepCellIdsByType(INTERP_KERNEL::NormalizedCellT * This method makes the assumption that da->getNumberOfTuples()getNumberOfCells(). This method makes the assumption that ids contained in 'da' * are in [0:getNumberOfCells()) */ -DataArrayInt *MEDCouplingUMesh::convertCellArrayPerGeoType(const DataArrayInt *da) const throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::convertCellArrayPerGeoType(const DataArrayInt *da) const { checkFullyDefined(); const int *conn=_nodal_connec->getConstPointer(); const int *connI=_nodal_connec_index->getConstPointer(); int nbOfCells=getNumberOfCells(); - std::set types=getAllTypes(); + std::set types(getAllGeoTypes()); int *tmp=new int[nbOfCells]; for(std::set::const_iterator iter=types.begin();iter!=types.end();iter++) { @@ -7002,7 +7160,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::keepSpecifiedCells(INTERP_KERNEL::Normalized * This method returns a vector of size 'this->getNumberOfCells()'. * This method retrieves for each cell in \a this if it is linear (false) or quadratic(true). */ -std::vector MEDCouplingUMesh::getQuadraticStatus() const throw(INTERP_KERNEL::Exception) +std::vector MEDCouplingUMesh::getQuadraticStatus() const { int ncell=getNumberOfCells(); std::vector ret(ncell); @@ -7071,7 +7229,7 @@ DataArrayDouble *MEDCouplingUMesh::getBarycenterAndOwner() const * \throw If \a this is not fully defined (coordinates and connectivity) * \throw If there is presence in nodal connectivity in \a this of node ids not in [0, \c this->getNumberOfNodes() ) */ -DataArrayDouble *MEDCouplingUMesh::computeIsoBarycenterOfNodesPerCell() const throw(INTERP_KERNEL::Exception) +DataArrayDouble *MEDCouplingUMesh::computeIsoBarycenterOfNodesPerCell() const { checkFullyDefined(); MEDCouplingAutoRefCountObjectPtr ret=DataArrayDouble::New(); @@ -7172,16 +7330,66 @@ DataArrayDouble *MEDCouplingUMesh::getPartBarycenterAndOwner(const int *begin, c return ret; } +/*! + * Returns a DataArrayDouble instance giving for each cell in \a this the equation of plane given by "a*X+b*Y+c*Z+d=0". + * So the returned instance will have 4 components and \c this->getNumberOfCells() tuples. + * So this method expects that \a this has a spaceDimension equal to 3 and meshDimension equal to 2. + * The computation of the plane equation is done using each time the 3 first nodes of 2D cells. + * This method is useful to detect 2D cells in 3D space that are not coplanar. + * + * \return DataArrayDouble * - a new instance of DataArrayDouble having 4 components and a number of tuples equal to number of cells in \a this. + * \throw If spaceDim!=3 or meshDim!=2. + * \throw If connectivity of \a this is invalid. + * \throw If connectivity of a cell in \a this points to an invalid node. + */ +DataArrayDouble *MEDCouplingUMesh::computePlaneEquationOf3DFaces() const +{ + MEDCouplingAutoRefCountObjectPtr ret(DataArrayDouble::New()); + int nbOfCells(getNumberOfCells()),nbOfNodes(getNumberOfNodes()); + if(getSpaceDimension()!=3 || getMeshDimension()!=2) + throw INTERP_KERNEL::Exception("MEDCouplingUMesh::computePlaneEquationOf3DFaces : This method must be applied on a mesh having meshDimension equal 2 and a spaceDimension equal to 3 !"); + ret->alloc(nbOfCells,4); + double *retPtr(ret->getPointer()); + const int *nodal(_nodal_connec->begin()),*nodalI(_nodal_connec_index->begin()); + const double *coor(_coords->begin()); + for(int i=0;i=3) + { + for(int j=0;j<3;j++) + { + int nodeId(nodal[nodalI[0]+1+j]); + if(nodeId>=0 && nodeIdcheckAllocated(); - MEDCouplingAutoRefCountObjectPtr ret=MEDCouplingUMesh::New(da->getName().c_str(),0); + MEDCouplingAutoRefCountObjectPtr ret=MEDCouplingUMesh::New(da->getName(),0); ret->setCoords(da); int nbOfTuples=da->getNumberOfTuples(); MEDCouplingAutoRefCountObjectPtr c=DataArrayInt::New(); @@ -7214,7 +7422,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::Build0DMeshFromCoords(DataArrayDouble *da) t * \throw If \a mesh1->getMeshDimension() < 0 or \a mesh2->getMeshDimension() < 0. * \throw If \a mesh1->getMeshDimension() != \a mesh2->getMeshDimension(). */ -MEDCouplingUMesh *MEDCouplingUMesh::MergeUMeshes(const MEDCouplingUMesh *mesh1, const MEDCouplingUMesh *mesh2) throw(INTERP_KERNEL::Exception) +MEDCouplingUMesh *MEDCouplingUMesh::MergeUMeshes(const MEDCouplingUMesh *mesh1, const MEDCouplingUMesh *mesh2) { std::vector tmp(2); tmp[0]=const_cast(mesh1); tmp[1]=const_cast(mesh2); @@ -7235,7 +7443,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::MergeUMeshes(const MEDCouplingUMesh *mesh1, * \throw If \a a[ *i* ]->getMeshDimension() < 0. * \throw If the meshes in \a a are of different dimension (getMeshDimension()). */ -MEDCouplingUMesh *MEDCouplingUMesh::MergeUMeshes(std::vector& a) throw(INTERP_KERNEL::Exception) +MEDCouplingUMesh *MEDCouplingUMesh::MergeUMeshes(std::vector& a) { std::size_t sz=a.size(); if(sz==0) @@ -7268,7 +7476,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::MergeUMeshes(std::vector& a) throw(INTERP_KERNEL::Exception) +MEDCouplingUMesh *MEDCouplingUMesh::MergeUMeshesLL(std::vector& a) { if(a.empty()) throw INTERP_KERNEL::Exception("MEDCouplingUMesh::MergeUMeshes : input array must be NON EMPTY !"); @@ -7340,7 +7548,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::MergeUMeshesLL(std::vectorgetMeshDimension() < 0 or \a mesh2->getMeshDimension() < 0. * \throw If \a mesh1->getMeshDimension() != \a mesh2->getMeshDimension(). */ -MEDCouplingUMesh *MEDCouplingUMesh::MergeUMeshesOnSameCoords(const MEDCouplingUMesh *mesh1, const MEDCouplingUMesh *mesh2) throw(INTERP_KERNEL::Exception) +MEDCouplingUMesh *MEDCouplingUMesh::MergeUMeshesOnSameCoords(const MEDCouplingUMesh *mesh1, const MEDCouplingUMesh *mesh2) { std::vector tmp(2); tmp[0]=mesh1; tmp[1]=mesh2; @@ -7454,7 +7662,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::FuseUMeshesOnSameCoords(const std::vectoralloc(curNbOfCells,1); std::copy(o2nPtr+offset,o2nPtr+offset+curNbOfCells,tmp->getPointer()); offset+=curNbOfCells; - tmp->setName(meshes[i]->getName().c_str()); + tmp->setName(meshes[i]->getName()); corr[i]=tmp; } return ret.retn(); @@ -7473,7 +7681,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::FuseUMeshesOnSameCoords(const std::vector& meshes) throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::PutUMeshesOnSameAggregatedCoords(const std::vector& meshes) { std::size_t sz=meshes.size(); if(sz==0 || sz==1) @@ -7527,7 +7735,7 @@ void MEDCouplingUMesh::PutUMeshesOnSameAggregatedCoords(const std::vector& meshes, double eps) throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::MergeNodesOnUMeshesSharingSameCoords(const std::vector& meshes, double eps) { if(meshes.empty()) return ; @@ -7755,7 +7963,7 @@ bool MEDCouplingUMesh::IsPyra5WellOriented(const int *begin, const int *end, con * \param [in] end end of nodal connectivity of a single polyhedron cell (excluded) * \param [out] res the result is put at the end of the vector without any alteration of the data. */ -void MEDCouplingUMesh::SimplifyPolyhedronCell(double eps, const DataArrayDouble *coords, const int *begin, const int *end, DataArrayInt *res) throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::SimplifyPolyhedronCell(double eps, const DataArrayDouble *coords, const int *begin, const int *end, DataArrayInt *res) { int nbFaces=std::count(begin+1,end,-1)+1; MEDCouplingAutoRefCountObjectPtr v=DataArrayDouble::New(); v->alloc(nbFaces,3); @@ -7847,7 +8055,7 @@ void MEDCouplingUMesh::SimplifyPolyhedronCell(double eps, const DataArrayDouble * \param [out] v the normalized vector of size 3 * \param [out] p the pos of plane */ -void MEDCouplingUMesh::ComputeVecAndPtOfFace(double eps, const double *coords, const int *begin, const int *end, double *v, double *p) throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::ComputeVecAndPtOfFace(double eps, const double *coords, const int *begin, const int *end, double *v, double *p) { std::size_t nbPoints=std::distance(begin,end); if(nbPoints<3) @@ -7893,7 +8101,7 @@ void MEDCouplingUMesh::ComputeVecAndPtOfFace(double eps, const double *coords, c * This method tries to obtain a well oriented polyhedron. * If the algorithm fails, an exception will be thrown. */ -void MEDCouplingUMesh::TryToCorrectPolyhedronOrientation(int *begin, int *end, const double *coords) throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::TryToCorrectPolyhedronOrientation(int *begin, int *end, const double *coords) { std::list< std::pair > edgesOK,edgesFinished; std::size_t nbOfFaces=std::count(begin,end,-1)+1; @@ -7970,7 +8178,7 @@ void MEDCouplingUMesh::TryToCorrectPolyhedronOrientation(int *begin, int *end, c * * \return a newly allocated array containing the connectivity of a polygon type enum included (NORM_POLYGON in pos#0) */ -DataArrayInt *MEDCouplingUMesh::buildUnionOf2DMesh() const throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::buildUnionOf2DMesh() const { if(getMeshDimension()!=2 || getSpaceDimension()!=2) throw INTERP_KERNEL::Exception("MEDCouplingUMesh::buildUnionOf2DMesh : meshdimension, spacedimension must be equal to 2 !"); @@ -8028,7 +8236,7 @@ DataArrayInt *MEDCouplingUMesh::buildUnionOf2DMesh() const throw(INTERP_KERNEL:: * * \return a newly allocated array containing the connectivity of a polygon type enum included (NORM_POLYHED in pos#0) */ -DataArrayInt *MEDCouplingUMesh::buildUnionOf3DMesh() const throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::buildUnionOf3DMesh() const { if(getMeshDimension()!=3 || getSpaceDimension()!=3) throw INTERP_KERNEL::Exception("MEDCouplingUMesh::buildUnionOf3DMesh : meshdimension, spacedimension must be equal to 2 !"); @@ -8053,7 +8261,7 @@ DataArrayInt *MEDCouplingUMesh::buildUnionOf3DMesh() const throw(INTERP_KERNEL:: * This method put in zip format into parameter 'zipFrmt' in full interlace mode. * This format is often asked by INTERP_KERNEL algorithms to avoid many indirections into coordinates array. */ -void MEDCouplingUMesh::FillInCompact3DMode(int spaceDim, int nbOfNodesInCell, const int *conn, const double *coo, double *zipFrmt) throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::FillInCompact3DMode(int spaceDim, int nbOfNodesInCell, const int *conn, const double *coo, double *zipFrmt) { double *w=zipFrmt; if(spaceDim==3) @@ -8071,12 +8279,12 @@ void MEDCouplingUMesh::FillInCompact3DMode(int spaceDim, int nbOfNodesInCell, co throw INTERP_KERNEL::Exception("MEDCouplingUMesh::FillInCompact3DMode : Invalid spaceDim specified : must be 2 or 3 !"); } -void MEDCouplingUMesh::writeVTKLL(std::ostream& ofs, const std::string& cellData, const std::string& pointData) const throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::writeVTKLL(std::ostream& ofs, const std::string& cellData, const std::string& pointData, DataArrayByte *byteData) const { int nbOfCells=getNumberOfCells(); if(nbOfCells<=0) throw INTERP_KERNEL::Exception("MEDCouplingUMesh::writeVTK : the unstructured mesh has no cells !"); - static const int PARAMEDMEM2VTKTYPETRADUCER[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,-1,4}; + static const int PARAMEDMEM2VTKTYPETRADUCER[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}; ofs << " <" << getVTKDataSetType() << ">\n"; ofs << " \n"; ofs << " \n" << pointData << std::endl; @@ -8085,11 +8293,11 @@ void MEDCouplingUMesh::writeVTKLL(std::ostream& ofs, const std::string& cellData ofs << " \n"; ofs << " \n"; if(getSpaceDimension()==3) - _coords->writeVTK(ofs,8,"Points"); + _coords->writeVTK(ofs,8,"Points",byteData); else { MEDCouplingAutoRefCountObjectPtr coo=_coords->changeNbOfComponents(3,0.); - coo->writeVTK(ofs,8,"Points"); + coo->writeVTK(ofs,8,"Points",byteData); } ofs << " \n"; ofs << " \n"; @@ -8120,12 +8328,12 @@ void MEDCouplingUMesh::writeVTKLL(std::ostream& ofs, const std::string& cellData } } types->transformWithIndArr(PARAMEDMEM2VTKTYPETRADUCER,PARAMEDMEM2VTKTYPETRADUCER+INTERP_KERNEL::NORM_MAXTYPE); - types->writeVTK(ofs,8,"UInt8","types"); - offsets->writeVTK(ofs,8,"Int32","offsets"); + types->writeVTK(ofs,8,"UInt8","types",byteData); + offsets->writeVTK(ofs,8,"Int32","offsets",byteData); if(szFaceOffsets!=0) {//presence of Polyhedra connectivity->reAlloc(szConn); - faceoffsets->writeVTK(ofs,8,"Int32","faceoffsets"); + faceoffsets->writeVTK(ofs,8,"Int32","faceoffsets",byteData); MEDCouplingAutoRefCountObjectPtr faces=DataArrayInt::New(); faces->alloc(szFaceOffsets,1); w1=faces->getPointer(); for(int i=0;iwriteVTK(ofs,8,"Int32","faces"); + faces->writeVTK(ofs,8,"Int32","faces",byteData); } - connectivity->writeVTK(ofs,8,"Int32","connectivity"); + connectivity->writeVTK(ofs,8,"Int32","connectivity",byteData); ofs << " \n"; ofs << " \n"; ofs << " \n"; } -void MEDCouplingUMesh::reprQuickOverview(std::ostream& stream) const throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::reprQuickOverview(std::ostream& stream) const { stream << "MEDCouplingUMesh C++ instance at " << this << ". Name : \"" << getName() << "\"."; if(_mesh_dim==-2) @@ -8175,14 +8383,17 @@ void MEDCouplingUMesh::reprQuickOverview(std::ostream& stream) const throw(INTER stream << std::endl << "Number of cells : " << lgth-1 << "."; } -std::string MEDCouplingUMesh::getVTKDataSetType() const throw(INTERP_KERNEL::Exception) +std::string MEDCouplingUMesh::getVTKDataSetType() const { return std::string("UnstructuredGrid"); } /*! * Partitions the first given 2D mesh using the second given 2D mesh as a tool, and - * returns a result mesh constituted by polygons. The meshes should be in 2D space. In + * returns a result mesh constituted by polygons. + * Thus the final result contains all nodes from m1 plus new nodes. However it doesn't necessarily contains + * all nodes from m2. + * The meshes should be in 2D space. In * addition, returns two arrays mapping cells of the result mesh to cells of the input * meshes. * \param [in] m1 - the first input mesh which is a partitioned object. @@ -8203,31 +8414,40 @@ std::string MEDCouplingUMesh::getVTKDataSetType() const throw(INTERP_KERNEL::Exc * \throw If the nodal connectivity of cells is not defined in any of the meshes. * \throw If any of the meshes is not a 2D mesh in 2D space. */ -MEDCouplingUMesh *MEDCouplingUMesh::Intersect2DMeshes(const MEDCouplingUMesh *m1, const MEDCouplingUMesh *m2, double eps, DataArrayInt *&cellNb1, DataArrayInt *&cellNb2) throw(INTERP_KERNEL::Exception) +MEDCouplingUMesh *MEDCouplingUMesh::Intersect2DMeshes(const MEDCouplingUMesh *m1, const MEDCouplingUMesh *m2, + double eps, DataArrayInt *&cellNb1, DataArrayInt *&cellNb2) { m1->checkFullyDefined(); m2->checkFullyDefined(); if(m1->getMeshDimension()!=2 || m1->getSpaceDimension()!=2 || m2->getMeshDimension()!=2 || m2->getSpaceDimension()!=2) throw INTERP_KERNEL::Exception("MEDCouplingUMesh::Intersect2DMeshes works on umeshes m1 AND m2 with meshdim equal to 2 and spaceDim equal to 2 too!"); + + // Step 1: compute all edge intersections (new nodes) std::vector< std::vector > intersectEdge1, colinear2, subDiv2; - MEDCouplingUMesh *m1Desc=0,*m2Desc=0; + MEDCouplingUMesh *m1Desc=0,*m2Desc=0; // descending connec. meshes DataArrayInt *desc1=0,*descIndx1=0,*revDesc1=0,*revDescIndx1=0,*desc2=0,*descIndx2=0,*revDesc2=0,*revDescIndx2=0; - std::vector addCoo,addCoordsQuadratic; + std::vector addCoo,addCoordsQuadratic; // coordinates of newly created nodes INTERP_KERNEL::QUADRATIC_PLANAR::_precision=eps; INTERP_KERNEL::QUADRATIC_PLANAR::_arc_detection_precision=eps; - IntersectDescending2DMeshes(m1,m2,eps,intersectEdge1,colinear2, subDiv2,m1Desc,desc1,descIndx1,revDesc1,revDescIndx1, - m2Desc,desc2,descIndx2,revDesc2,revDescIndx2,addCoo); + IntersectDescending2DMeshes(m1,m2,eps,intersectEdge1,colinear2, subDiv2, + m1Desc,desc1,descIndx1,revDesc1,revDescIndx1, + addCoo, m2Desc,desc2,descIndx2,revDesc2,revDescIndx2); revDesc1->decrRef(); revDescIndx1->decrRef(); revDesc2->decrRef(); revDescIndx2->decrRef(); MEDCouplingAutoRefCountObjectPtr dd1(desc1),dd2(descIndx1),dd3(desc2),dd4(descIndx2); MEDCouplingAutoRefCountObjectPtr dd5(m1Desc),dd6(m2Desc); + + // Step 2: re-order newly created nodes according to the ordering found in m2 std::vector< std::vector > intersectEdge2; BuildIntersectEdges(m1Desc,m2Desc,addCoo,subDiv2,intersectEdge2); subDiv2.clear(); dd5=0; dd6=0; + + // Step 3: std::vector cr,crI; //no DataArrayInt because interface with Geometric2D std::vector cNb1,cNb2; //no DataArrayInt because interface with Geometric2D BuildIntersecting2DCellsFromEdges(eps,m1,desc1->getConstPointer(),descIndx1->getConstPointer(),intersectEdge1,colinear2,m2,desc2->getConstPointer(),descIndx2->getConstPointer(),intersectEdge2,addCoo, /* outputs -> */addCoordsQuadratic,cr,crI,cNb1,cNb2); - // + + // Step 4: Prepare final result: MEDCouplingAutoRefCountObjectPtr addCooDa=DataArrayDouble::New(); addCooDa->alloc((int)(addCoo.size())/2,2); std::copy(addCoo.begin(),addCoo.end(),addCooDa->getPointer()); @@ -8248,6 +8468,15 @@ MEDCouplingUMesh *MEDCouplingUMesh::Intersect2DMeshes(const MEDCouplingUMesh *m1 return ret.retn(); } + +/** + * Private. Third step of the partitioning algorithm (Intersect2DMeshes): reconstruct full 2D cells from the + * (newly created) nodes corresponding to the edge intersections. + * Output params: + * @param[out] cr, crI connectivity of the resulting mesh + * @param[out] cNb1, cNb2 correspondance arrays giving for the merged mesh the initial cells IDs in m1 / m2 + * TODO: describe input parameters + */ void MEDCouplingUMesh::BuildIntersecting2DCellsFromEdges(double eps, const MEDCouplingUMesh *m1, const int *desc1, const int *descIndx1, const std::vector >& intesctEdges1, const std::vector< std::vector >& colinear2, const MEDCouplingUMesh *m2, const int *desc2, const int *descIndx2, const std::vector >& intesctEdges2, @@ -8266,6 +8495,7 @@ void MEDCouplingUMesh::BuildIntersecting2DCellsFromEdges(double eps, const MEDCo int offset3=offset2+((int)addCoords.size())/2; MEDCouplingAutoRefCountObjectPtr bbox1Arr(m1->getBoundingBoxForBBTree()),bbox2Arr(m2->getBoundingBoxForBBTree()); const double *bbox1(bbox1Arr->begin()),*bbox2(bbox2Arr->begin()); + // Here a BBTree on 2D-cells, not on segments: BBTree myTree(bbox2,0,0,m2->getNumberOfCells(),eps); int ncell1=m1->getNumberOfCells(); crI.push_back(0); @@ -8278,7 +8508,9 @@ void MEDCouplingUMesh::BuildIntersecting2DCellsFromEdges(double eps, const MEDCo INTERP_KERNEL::QuadraticPolygon pol1; INTERP_KERNEL::NormalizedCellType typ=(INTERP_KERNEL::NormalizedCellType)conn1[connI1[i]]; const INTERP_KERNEL::CellModel& cm=INTERP_KERNEL::CellModel::GetCellModel(typ); + // Populate mapp and mappRev with nodes from the current cell (i) from mesh1 - this also builds the Node* objects: MEDCouplingUMeshBuildQPFromMesh3(coo1,offset1,coo2,offset2,addCoords,desc1+descIndx1[i],desc1+descIndx1[i+1],intesctEdges1,/* output */mapp,mappRev); + // pol1 is the full cell from mesh2, in QP format, with all the additional intersecting nodes. pol1.buildFromCrudeDataArray(mappRev,cm.isQuadratic(),conn1+connI1[i]+1,coo1, desc1+descIndx1[i],desc1+descIndx1[i+1],intesctEdges1); // @@ -8288,16 +8520,18 @@ void MEDCouplingUMesh::BuildIntersecting2DCellsFromEdges(double eps, const MEDCo for(it1.first();!it1.finished();it1.next()) edges1.insert(it1.current()->getPtr()); // - std::map > edgesIn2ForShare; + std::map > edgesIn2ForShare; // common edges std::vector pol2s(candidates2.size()); int ii=0; for(std::vector::const_iterator it2=candidates2.begin();it2!=candidates2.end();it2++,ii++) { INTERP_KERNEL::NormalizedCellType typ2=(INTERP_KERNEL::NormalizedCellType)conn2[connI2[*it2]]; const INTERP_KERNEL::CellModel& cm2=INTERP_KERNEL::CellModel::GetCellModel(typ2); + // Complete mapping with elements coming from the current cell it2 in mesh2: MEDCouplingUMeshBuildQPFromMesh3(coo1,offset1,coo2,offset2,addCoords,desc2+descIndx2[*it2],desc2+descIndx2[*it2+1],intesctEdges2,/* output */mapp,mappRev); + // pol2 is the new QP in the final merged result. pol2s[ii].buildFromCrudeDataArray2(mappRev,cm2.isQuadratic(),conn2+connI2[*it2]+1,coo2,desc2+descIndx2[*it2],desc2+descIndx2[*it2+1],intesctEdges2, - pol1,desc1+descIndx1[i],desc1+descIndx1[i+1],intesctEdges1,colinear2,edgesIn2ForShare); + pol1,desc1+descIndx1[i],desc1+descIndx1[i+1],intesctEdges1,colinear2, /* output */ edgesIn2ForShare); } ii=0; for(std::vector::const_iterator it2=candidates2.begin();it2!=candidates2.end();it2++,ii++) @@ -8307,6 +8541,8 @@ void MEDCouplingUMesh::BuildIntersecting2DCellsFromEdges(double eps, const MEDCo //MEDCouplingUMeshAssignOnLoc(pol1,pol2,desc1+descIndx1[i],desc1+descIndx1[i+1],intesctEdges1,desc2+descIndx2[*it2],desc2+descIndx2[*it2+1],intesctEdges2,colinear2); pol1.buildPartitionsAbs(pol2s[ii],edges1,edgesBoundary2,mapp,i,*it2,offset3,addCoordsQuadratic,cr,crI,cNb1,cNb2); } + // Deals with remaining (non-consumed) edges from m1: these are the edges that were never touched + // by m2 but that we still want to keep in the final result. if(!edges1.empty()) { try @@ -8326,15 +8562,19 @@ void MEDCouplingUMesh::BuildIntersecting2DCellsFromEdges(double eps, const MEDCo /*! * This method is private and is the first step of Partition of 2D mesh (spaceDim==2 and meshDim==2). - * + * It builds the descending connectivity of the two meshes, and then using a binary tree + * it computes the edge intersections. This results in new points being created : they're stored in addCoo. + * Documentation about parameters colinear2 and subDiv2 can be found in method QuadraticPolygon::splitAbs(). */ void MEDCouplingUMesh::IntersectDescending2DMeshes(const MEDCouplingUMesh *m1, const MEDCouplingUMesh *m2, double eps, std::vector< std::vector >& intersectEdge1, std::vector< std::vector >& colinear2, std::vector< std::vector >& subDiv2, MEDCouplingUMesh *& m1Desc, DataArrayInt *&desc1, DataArrayInt *&descIndx1, DataArrayInt *&revDesc1, DataArrayInt *&revDescIndx1, - MEDCouplingUMesh *& m2Desc, DataArrayInt *&desc2, DataArrayInt *&descIndx2, DataArrayInt *&revDesc2, DataArrayInt *&revDescIndx2, - std::vector& addCoo) throw(INTERP_KERNEL::Exception) + std::vector& addCoo, + MEDCouplingUMesh *& m2Desc, DataArrayInt *&desc2, DataArrayInt *&descIndx2, DataArrayInt *&revDesc2, DataArrayInt *&revDescIndx2) + throw(INTERP_KERNEL::Exception) { static const int SPACEDIM=2; + // Build desc connectivity desc1=DataArrayInt::New(); descIndx1=DataArrayInt::New(); revDesc1=DataArrayInt::New(); revDescIndx1=DataArrayInt::New(); desc2=DataArrayInt::New(); descIndx2=DataArrayInt::New(); @@ -8347,27 +8587,42 @@ void MEDCouplingUMesh::IntersectDescending2DMeshes(const MEDCouplingUMesh *m1, c MEDCouplingAutoRefCountObjectPtr dd9(m1Desc),dd10(m2Desc); const int *c1=m1Desc->getNodalConnectivity()->getConstPointer(); const int *ci1=m1Desc->getNodalConnectivityIndex()->getConstPointer(); + + // Build BB tree of all edges in the tool mesh (second mesh) MEDCouplingAutoRefCountObjectPtr bbox1Arr(m1Desc->getBoundingBoxForBBTree()),bbox2Arr(m2Desc->getBoundingBoxForBBTree()); const double *bbox1(bbox1Arr->begin()),*bbox2(bbox2Arr->begin()); - int ncell1=m1Desc->getNumberOfCells(); - int ncell2=m2Desc->getNumberOfCells(); - intersectEdge1.resize(ncell1); - colinear2.resize(ncell2); - subDiv2.resize(ncell2); + int nDescCell1=m1Desc->getNumberOfCells(); + int nDescCell2=m2Desc->getNumberOfCells(); + intersectEdge1.resize(nDescCell1); + colinear2.resize(nDescCell2); + subDiv2.resize(nDescCell2); BBTree myTree(bbox2,0,0,m2Desc->getNumberOfCells(),-eps); + std::vector candidates1(1); int offset1=m1->getNumberOfNodes(); int offset2=offset1+m2->getNumberOfNodes(); - for(int i=0;i candidates2; + std::vector candidates2; // edges of mesh2 candidate for intersection myTree.getIntersectingElems(bbox1+i*2*SPACEDIM,candidates2); - if(!candidates2.empty()) + if(!candidates2.empty()) // candidates2 holds edges from the second mesh potentially intersecting current edge i in mesh1 { std::map map1,map2; + // pol2 is not necessarily a closed polygon: just a set of (quadratic) edges (same as candidates2) in the Geometric DS format INTERP_KERNEL::QuadraticPolygon *pol2=MEDCouplingUMeshBuildQPFromMesh(m2Desc,candidates2,map2); candidates1[0]=i; INTERP_KERNEL::QuadraticPolygon *pol1=MEDCouplingUMeshBuildQPFromMesh(m1Desc,candidates1,map1); + // This following part is to avoid that some removed nodes (for example due to a merge between pol1 and pol2) are replaced by a newly created one + // This trick guarantees that Node * are discriminant (i.e. form a unique identifier) + std::set nodes; + pol1->getAllNodes(nodes); pol2->getAllNodes(nodes); + std::size_t szz(nodes.size()); + std::vector< MEDCouplingAutoRefCountObjectPtr > nodesSafe(szz); + std::set::const_iterator itt(nodes.begin()); + for(std::size_t iii=0;iiiincrRef(); nodesSafe[iii]=*itt; } + // end of protection + // Performs egde cutting: pol1->splitAbs(*pol2,map1,map2,offset1,offset2,candidates2,intersectEdge1[i],i,colinear2,subDiv2,addCoo); delete pol2; delete pol1; @@ -8384,13 +8639,18 @@ void MEDCouplingUMesh::IntersectDescending2DMeshes(const MEDCouplingUMesh *m1, c * This method has 4 inputs : * - a mesh 'm1' with meshDim==1 and a SpaceDim==2 * - a mesh 'm2' with meshDim==1 and a SpaceDim==2 - * - subDiv of size 'm2->getNumberOfCells()' that lists for each seg cell in 'm' the splitting node ids in randomly sorted. - * The aim of this method is to sort the splitting nodes, if any, and to put in 'intersectEdge' output paramter based on edges of mesh 'm2' - * \param m1 is expected to be a mesh of meshDimension equal to 1 and spaceDim equal to 2. No check of that is performed by this method. Only present for its coords in case of 'subDiv' shares some nodes of 'm1' + * - subDiv of size 'm2->getNumberOfCells()' that lists for each seg cell in 'm' the splitting node ids randomly sorted. + * The aim of this method is to sort the splitting nodes, if any, and to put them in 'intersectEdge' output parameter based on edges of mesh 'm2' + * Nodes end up lying consecutively on a cutted edge. + * \param m1 is expected to be a mesh of meshDimension equal to 1 and spaceDim equal to 2. No check of that is performed by this method. + * (Only present for its coords in case of 'subDiv' shares some nodes of 'm1') * \param m2 is expected to be a mesh of meshDimension equal to 1 and spaceDim equal to 2. No check of that is performed by this method. - * \param addCoo input parameter with additionnal nodes linked to intersection of the 2 meshes. + * \param addCoo input parameter with additional nodes linked to intersection of the 2 meshes. + * \param[out] intersectEdge the same content as subDiv, but correclty oriented. */ -void MEDCouplingUMesh::BuildIntersectEdges(const MEDCouplingUMesh *m1, const MEDCouplingUMesh *m2, const std::vector& addCoo, const std::vector< std::vector >& subDiv, std::vector< std::vector >& intersectEdge) throw(INTERP_KERNEL::Exception) +void MEDCouplingUMesh::BuildIntersectEdges(const MEDCouplingUMesh *m1, const MEDCouplingUMesh *m2, + const std::vector& addCoo, + const std::vector< std::vector >& subDiv, std::vector< std::vector >& intersectEdge) { int offset1=m1->getNumberOfNodes(); int ncell=m2->getNumberOfCells(); @@ -8607,7 +8867,7 @@ void MEDCouplingUMesh::assemblyForSplitFrom3DSurf(const std::vector< std::pair=4) @@ -8706,7 +8966,7 @@ bool MEDCouplingUMesh::BuildConvexEnvelopOf2DCellJarvis(const double *coords, co * \param [in] offsetForRemoval (by default 0) offset so that for each i in [0,arrIndx->getNumberOfTuples()-1) removal process will be performed in the following range [arr+arrIndx[i]+offsetForRemoval,arr+arr[i+1]) * \return true if \b arr and \b arrIndx have been modified, false if not. */ -bool MEDCouplingUMesh::RemoveIdsFromIndexedArrays(const int *idsToRemoveBg, const int *idsToRemoveEnd, DataArrayInt *arr, DataArrayInt *arrIndx, int offsetForRemoval) throw(INTERP_KERNEL::Exception) +bool MEDCouplingUMesh::RemoveIdsFromIndexedArrays(const int *idsToRemoveBg, const int *idsToRemoveEnd, DataArrayInt *arr, DataArrayInt *arrIndx, int offsetForRemoval) { if(!arrIndx || !arr) throw INTERP_KERNEL::Exception("MEDCouplingUMesh::RemoveIdsFromIndexedArrays : some input arrays are empty !"); @@ -8847,7 +9107,7 @@ void MEDCouplingUMesh::ExtractFromIndexedArrays2(int idsOfSelectStart, int idsOf int *work=arrIo->getPointer(); *work++=0; int lgth=0; - for(std::size_t i=0;i=0 && idsItalloc(lgth,1); work=arro->getPointer(); idsIt=idsOfSelectStart; - for(std::size_t i=0;i=0 && arrIndxPtr[idsIt+1]<=maxSizeOfArr) work=std::copy(arrInPtr+arrIndxPtr[idsIt],arrInPtr+arrIndxPtr[idsIt+1],work); @@ -9007,7 +9267,7 @@ void MEDCouplingUMesh::SetPartOfIndexedArraysSameIdx(const int *idsOfSelectBg, c * \return a newly allocated DataArray that stores all ids fetched by the gradually spread process. * \sa MEDCouplingUMesh::ComputeSpreadZoneGraduallyFromSeed, MEDCouplingUMesh::partitionBySpreadZone */ -DataArrayInt *MEDCouplingUMesh::ComputeSpreadZoneGradually(const DataArrayInt *arrIn, const DataArrayInt *arrIndxIn) throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::ComputeSpreadZoneGradually(const DataArrayInt *arrIn, const DataArrayInt *arrIndxIn) { int seed=0,nbOfDepthPeelingPerformed=0; return ComputeSpreadZoneGraduallyFromSeed(&seed,&seed+1,arrIn,arrIndxIn,-1,nbOfDepthPeelingPerformed); @@ -9029,7 +9289,7 @@ DataArrayInt *MEDCouplingUMesh::ComputeSpreadZoneGradually(const DataArrayInt *a * \return a newly allocated DataArray that stores all ids fetched by the gradually spread process. * \sa MEDCouplingUMesh::partitionBySpreadZone */ -DataArrayInt *MEDCouplingUMesh::ComputeSpreadZoneGraduallyFromSeed(const int *seedBg, const int *seedEnd, const DataArrayInt *arrIn, const DataArrayInt *arrIndxIn, int nbOfDepthPeeling, int& nbOfDepthPeelingPerformed) throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::ComputeSpreadZoneGraduallyFromSeed(const int *seedBg, const int *seedEnd, const DataArrayInt *arrIn, const DataArrayInt *arrIndxIn, int nbOfDepthPeeling, int& nbOfDepthPeelingPerformed) { nbOfDepthPeelingPerformed=0; if(!arrIndxIn) @@ -9045,7 +9305,7 @@ DataArrayInt *MEDCouplingUMesh::ComputeSpreadZoneGraduallyFromSeed(const int *se return ComputeSpreadZoneGraduallyFromSeedAlg(fetched,seedBg,seedEnd,arrIn,arrIndxIn,nbOfDepthPeeling,nbOfDepthPeelingPerformed); } -DataArrayInt *MEDCouplingUMesh::ComputeSpreadZoneGraduallyFromSeedAlg(std::vector& fetched, const int *seedBg, const int *seedEnd, const DataArrayInt *arrIn, const DataArrayInt *arrIndxIn, int nbOfDepthPeeling, int& nbOfDepthPeelingPerformed) throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::ComputeSpreadZoneGraduallyFromSeedAlg(std::vector& fetched, const int *seedBg, const int *seedEnd, const DataArrayInt *arrIn, const DataArrayInt *arrIndxIn, int nbOfDepthPeeling, int& nbOfDepthPeelingPerformed) { nbOfDepthPeelingPerformed=0; if(!seedBg || !seedEnd || !arrIn || !arrIndxIn) @@ -9209,7 +9469,7 @@ void MEDCouplingUMesh::SetPartOfIndexedArraysSameIdx2(int start, int end, int st * * \return a newly allocated mesh lying on the same coords than \b this with same meshdimension than \b this. */ -MEDCouplingUMesh *MEDCouplingUMesh::buildSpreadZonesWithPoly() const throw(INTERP_KERNEL::Exception) +MEDCouplingUMesh *MEDCouplingUMesh::buildSpreadZonesWithPoly() const { checkFullyDefined(); int mdim=getMeshDimension(); @@ -9219,7 +9479,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::buildSpreadZonesWithPoly() const throw(INTER std::vector partition=partitionBySpreadZone(); std::vector< MEDCouplingAutoRefCountObjectPtr > partitionAuto; partitionAuto.reserve(partition.size()); std::copy(partition.begin(),partition.end(),std::back_insert_iterator > >(partitionAuto)); - MEDCouplingAutoRefCountObjectPtr ret=MEDCouplingUMesh::New(getName().c_str(),mdim); + MEDCouplingAutoRefCountObjectPtr ret=MEDCouplingUMesh::New(getName(),mdim); ret->setCoords(getCoords()); ret->allocateCells((int)partition.size()); // @@ -9251,7 +9511,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::buildSpreadZonesWithPoly() const throw(INTER * This method only needs a well defined connectivity. Coordinates are not considered here. * This method returns a vector of \b newly allocated arrays that the caller has to deal with. */ -std::vector MEDCouplingUMesh::partitionBySpreadZone() const throw(INTERP_KERNEL::Exception) +std::vector MEDCouplingUMesh::partitionBySpreadZone() const { int nbOfCellsCur=getNumberOfCells(); std::vector ret; @@ -9282,7 +9542,7 @@ std::vector MEDCouplingUMesh::partitionBySpreadZone() const thro * \return a newly allocated DataArrayInt to be managed by the caller. * \throw In case of \a code has not the right format (typically of size 3*n) */ -DataArrayInt *MEDCouplingUMesh::ComputeRangesFromTypeDistribution(const std::vector& code) throw(INTERP_KERNEL::Exception) +DataArrayInt *MEDCouplingUMesh::ComputeRangesFromTypeDistribution(const std::vector& code) { MEDCouplingAutoRefCountObjectPtr ret=DataArrayInt::New(); std::size_t nb=code.size()/3; @@ -9314,19 +9574,22 @@ DataArrayInt *MEDCouplingUMesh::ComputeRangesFromTypeDistribution(const std::vec * 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. + * + * \warning This method operates on each cells in this independantly ! So it can leads to non conform mesh in returned value ! If you expect to have a conform mesh in output + * the policy PLANAR_FACE_6 should be used on a mesh sorted with MEDCoupling1SGTUMesh::sortHexa8EachOther. * * \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 + * \sa MEDCouplingUMesh::simplexize, MEDCoupling1SGTUMesh::sortHexa8EachOther */ -MEDCoupling1SGTUMesh *MEDCouplingUMesh::tetrahedrize(int policy, DataArrayInt *& n2oCells, int& nbOfAdditionalPoints) const throw(INTERP_KERNEL::Exception) +MEDCoupling1SGTUMesh *MEDCouplingUMesh::tetrahedrize(int policy, DataArrayInt *& n2oCells, int& nbOfAdditionalPoints) const { 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 ret0(MEDCoupling1SGTUMesh::New(getName().c_str(),INTERP_KERNEL::NORM_TETRA4)); + MEDCouplingAutoRefCountObjectPtr ret0(MEDCoupling1SGTUMesh::New(getName(),INTERP_KERNEL::NORM_TETRA4)); MEDCouplingAutoRefCountObjectPtr ret(DataArrayInt::New()); ret->alloc(nbOfCells,1); int *retPt(ret->getPointer()); MEDCouplingAutoRefCountObjectPtr newConn(DataArrayInt::New()); newConn->alloc(0,1);