#include "BBTreeDst.txx"
#include "SplitterTetra.hxx"
#include "DirectedBoundingBox.hxx"
+#include "InterpKernelMatrixTools.hxx"
#include "InterpKernelMeshQuality.hxx"
#include "InterpKernelCellSimplify.hxx"
#include "InterpKernelGeo2DEdgeArcCircle.hxx"
* 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)
{
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<int> 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<int> 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<int> 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);
+ }
}
}
}
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<DataArrayDouble> 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<nbOfCells;i++,nodalI++,retPtr+=4)
+ {
+ double matrix[16]={0,0,0,1,0,0,0,1,0,0,0,1,1,1,1,0},matrix2[16];
+ if(nodalI[1]-nodalI[0]>=3)
+ {
+ for(int j=0;j<3;j++)
+ {
+ int nodeId(nodal[nodalI[0]+1+j]);
+ if(nodeId>=0 && nodeId<nbOfNodes)
+ std::copy(coor+nodeId*3,coor+(nodeId+1)*3,matrix+4*j);
+ else
+ {
+ std::ostringstream oss; oss << "MEDCouplingUMesh::computePlaneEquationOf3DFaces : invalid 2D cell #" << i << " ! This cell points to an invalid nodeId : " << nodeId << " !";
+ throw INTERP_KERNEL::Exception(oss.str().c_str());
+ }
+ }
+ }
+ else
+ {
+ std::ostringstream oss; oss << "MEDCouplingUMesh::computePlaneEquationOf3DFaces : invalid 2D cell #" << i << " ! Must be constitued by more than 3 nodes !";
+ throw INTERP_KERNEL::Exception(oss.str().c_str());
+ }
+ INTERP_KERNEL::inverseMatrix(matrix,4,matrix2);
+ retPtr[0]=matrix2[3]; retPtr[1]=matrix2[7]; retPtr[2]=matrix2[11]; retPtr[3]=matrix2[15];
+ }
+ return ret.retn();
+}
+
/*!
* This method expects as input a DataArrayDouble non nul instance 'da' that should be allocated. If not an exception is thrown.
*
MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> dd5(m1Desc),dd6(m2Desc);
std::vector< std::vector<int> > intersectEdge2;
BuildIntersectEdges(m1Desc,m2Desc,addCoo,subDiv2,intersectEdge2);
- subDiv2.clear(); //dd5=0; dd6=0;
+ subDiv2.clear(); dd5=0; dd6=0;
std::vector<int> cr,crI; //no DataArrayInt because interface with Geometric2D
std::vector<int> cNb1,cNb2; //no DataArrayInt because interface with Geometric2D
BuildIntersecting2DCellsFromEdges(eps,m1,desc1->getConstPointer(),descIndx1->getConstPointer(),intersectEdge1,colinear2,m2,desc2->getConstPointer(),descIndx2->getConstPointer(),intersectEdge2,addCoo,
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 a some remove nodes (for example due to a merge between pol1 and pol2) can be replaced by a newlt created one
+ // This trick garanties that Node * are discriminant
+ std::set<INTERP_KERNEL::Node *> nodes;
+ pol1->getAllNodes(nodes); pol2->getAllNodes(nodes);
+ std::size_t szz(nodes.size());
+ std::vector< MEDCouplingAutoRefCountObjectPtr<INTERP_KERNEL::Node> > nodesSafe(szz);
+ std::set<INTERP_KERNEL::Node *>::const_iterator itt(nodes.begin());
+ for(std::size_t iii=0;iii<szz;iii++,itt++)
+ { (*itt)->incrRef(); nodesSafe[iii]=*itt; }
+ // end of protection
pol1->splitAbs(*pol2,map1,map2,offset1,offset2,candidates2,intersectEdge1[i],i,colinear2,subDiv2,addCoo);
delete pol2;
delete pol1;
* 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
{
