From b6e0752eacaf3a3905d3a9013417802b175d27e7 Mon Sep 17 00:00:00 2001 From: jfa Date: Fri, 19 Jul 2024 12:53:33 +0100 Subject: [PATCH] Avoid anonymous namespaces and 'using' instructions in headers. Simplify tests code. Continue refactoring of Hexahedron::compute() method. --- .../HexahedronCanonicalShapesTest.cxx | 69 +- src/StdMeshers.test/HexahedronTest.cxx | 120 ++- src/StdMeshers/StdMeshers_Cartesian_3D.cxx | 7 +- .../StdMeshers_Cartesian_3D_Grid.cxx | 359 +++++++- .../StdMeshers_Cartesian_3D_Grid.hxx | 489 ++--------- .../StdMeshers_Cartesian_3D_Hexahedron.cxx | 229 ++--- .../StdMeshers_Cartesian_3D_Hexahedron.hxx | 790 +++++++++--------- test/data/HexahedronTest/NRTMJ4.brep | 414 +++------ 8 files changed, 1137 insertions(+), 1340 deletions(-) diff --git a/src/StdMeshers.test/HexahedronCanonicalShapesTest.cxx b/src/StdMeshers.test/HexahedronCanonicalShapesTest.cxx index fe3ab4c95..3f4b16f56 100644 --- a/src/StdMeshers.test/HexahedronCanonicalShapesTest.cxx +++ b/src/StdMeshers.test/HexahedronCanonicalShapesTest.cxx @@ -39,11 +39,12 @@ #include #include +using namespace StdMeshers::Cartesian3D; + // Helper functions! // Build Grid // Require building mesh -// Require building shape. For test load shapes from memory in .brep files seems the simplest -// +// Require building shape. /*! * \brief Mock mesh @@ -69,62 +70,27 @@ struct CartesianHypo: public StdMeshers_CartesianParameters3D /*! * \brief Initialize the grid and intesersectors of grid with the geometry */ -void GridInitAndIntersectWithShape( Grid& grid, +void GridInitAndIntersectWithShape (Grid& grid, double gridSpacing, double theSizeThreshold, const TopoDS_Shape theShape, - std::map< TGeomID, vector< TGeomID > >& edge2faceIDsMap, - const int /*numOfThreads*/ ) + TEdge2faceIDsMap& edge2faceIDsMap, + const int theNumOfThreads) { - std::vector< TopoDS_Shape > faceVec; - TopTools_MapOfShape faceMap; - TopExp_Explorer fExp; - for ( fExp.Init( theShape, TopAbs_FACE ); fExp.More(); fExp.Next() ) - { - bool isNewFace = faceMap.Add( fExp.Current() ); - if ( !grid._toConsiderInternalFaces ) - if ( !isNewFace || fExp.Current().Orientation() == TopAbs_INTERNAL ) - // remove an internal face - faceMap.Remove( fExp.Current() ); - } - faceVec.reserve( faceMap.Extent() ); - faceVec.assign( faceMap.cbegin(), faceMap.cend() ); - - vector facesItersectors( faceVec.size() ); - - Bnd_Box shapeBox; - for ( size_t i = 0; i < faceVec.size(); ++i ) - { - facesItersectors[i]._face = TopoDS::Face( faceVec[i] ); - facesItersectors[i]._faceID = grid.ShapeID( faceVec[i] ); - facesItersectors[i]._grid = &grid; - shapeBox.Add( facesItersectors[i].GetFaceBndBox() ); - } // Canonical axes(i,j,k) double axisDirs[9] = {1.,0.,0.,0.,1.,0.,0.,0.,1.}; - - Tools::GetExactBndBox( faceVec, axisDirs, shapeBox ); - vector xCoords, yCoords, zCoords; - std::unique_ptr myHypo( new CartesianHypo() ); std::vector grdSpace = { std::to_string(gridSpacing) }; std::vector intPnts; - myHypo->SetGridSpacing(grdSpace, intPnts, 0 ); // Spacing in dir 0 - myHypo->SetGridSpacing(grdSpace, intPnts, 1 ); // Spacing in dir 1 - myHypo->SetGridSpacing(grdSpace, intPnts, 2 ); // Spacing in dir 2 - myHypo->SetSizeThreshold(theSizeThreshold); // set threshold - myHypo->GetCoordinates(xCoords, yCoords, zCoords, shapeBox); - grid.SetCoordinates( xCoords, yCoords, zCoords, axisDirs, shapeBox ); - - for ( size_t i = 0; i < facesItersectors.size(); ++i ) - facesItersectors[i].Intersect(); - - for ( size_t i = 0; i < facesItersectors.size(); ++i ) - facesItersectors[i].StoreIntersections(); - - grid.ComputeNodes( *grid._helper ); - grid.GetEdgesToImplement( edge2faceIDsMap, theShape, faceVec ); -} + std::unique_ptr aHypo ( new CartesianHypo() ); + aHypo->SetGridSpacing(grdSpace, intPnts, 0 ); // Spacing in dir 0 + aHypo->SetGridSpacing(grdSpace, intPnts, 1 ); // Spacing in dir 1 + aHypo->SetGridSpacing(grdSpace, intPnts, 2 ); // Spacing in dir 2 + aHypo->SetSizeThreshold(theSizeThreshold); // set threshold + aHypo->SetAxisDirs(axisDirs); + + grid.GridInitAndInterserctWithShape(theShape, edge2faceIDsMap, aHypo.get(), theNumOfThreads, false); +} /*! * \brief Test runner @@ -139,6 +105,7 @@ bool testShape (const TopoDS_Shape theShape, std::unique_ptr aMesh( new SMESH_Mesh_Test() ); aMesh->ShapeToMesh( theShape ); SMESH_MesherHelper helper( *aMesh ); + Grid grid; grid._helper = &helper; grid._toAddEdges = toAddEdges; @@ -147,11 +114,11 @@ bool testShape (const TopoDS_Shape theShape, grid._toUseThresholdForInternalFaces = false; grid._toUseQuanta = false; grid._sizeThreshold = theSizeThreshold; - grid.InitGeometry( theShape ); - std::map< TGeomID, vector< TGeomID > > edge2faceIDsMap; + TEdge2faceIDsMap edge2faceIDsMap; GridInitAndIntersectWithShape( grid, theGridSpacing, theSizeThreshold, theShape, edge2faceIDsMap, 1 ); + Hexahedron hex( &grid ); int nbAdded = hex.MakeElements( helper, edge2faceIDsMap, 1 ); if (nbAdded != theNbCreatedExpected) { diff --git a/src/StdMeshers.test/HexahedronTest.cxx b/src/StdMeshers.test/HexahedronTest.cxx index e993a96f5..dd78c796e 100644 --- a/src/StdMeshers.test/HexahedronTest.cxx +++ b/src/StdMeshers.test/HexahedronTest.cxx @@ -39,11 +39,12 @@ #include #include +using namespace StdMeshers::Cartesian3D; + // Helper functions! // Build Grid // Require building mesh -// Require building shape. For test load shapes from memory in .brep files seems the simplest -// +// Require building shape. /*! * \brief Mock mesh @@ -78,60 +79,26 @@ void loadBrepShape( std::string shapeName, TopoDS_Shape & shape ) /*! * \brief Initialize the grid and intesersectors of grid with the geometry */ -void GridInitAndIntersectWithShape( Grid& grid, +void GridInitAndIntersectWithShape (Grid& grid, double gridSpacing, double theSizeThreshold, const TopoDS_Shape theShape, - std::map< TGeomID, vector< TGeomID > >& edge2faceIDsMap, - const int /*numOfThreads*/ ) + TEdge2faceIDsMap& edge2faceIDsMap, + const int theNumOfThreads) { - std::vector< TopoDS_Shape > faceVec; - TopTools_MapOfShape faceMap; - TopExp_Explorer fExp; - for ( fExp.Init( theShape, TopAbs_FACE ); fExp.More(); fExp.Next() ) - { - bool isNewFace = faceMap.Add( fExp.Current() ); - if ( !grid._toConsiderInternalFaces ) - if ( !isNewFace || fExp.Current().Orientation() == TopAbs_INTERNAL ) - // remove an internal face - faceMap.Remove( fExp.Current() ); - } - faceVec.reserve( faceMap.Extent() ); - faceVec.assign( faceMap.cbegin(), faceMap.cend() ); - - vector facesItersectors( faceVec.size() ); - - Bnd_Box shapeBox; - for ( size_t i = 0; i < faceVec.size(); ++i ) - { - facesItersectors[i]._face = TopoDS::Face( faceVec[i] ); - facesItersectors[i]._faceID = grid.ShapeID( faceVec[i] ); - facesItersectors[i]._grid = &grid; - shapeBox.Add( facesItersectors[i].GetFaceBndBox() ); - } // Canonical axes(i,j,k) - double axisDirs[9] = {1.,0.,0.,0.,1.,0.,0.,0.,1.}; - - Tools::GetExactBndBox( faceVec, axisDirs, shapeBox ); - vector xCoords, yCoords, zCoords; - std::unique_ptr myHypo( new CartesianHypo() ); + double axisDirs[9] = {1.,0.,0., 0.,1.,0., 0.,0.,1.}; std::vector grdSpace = { std::to_string(gridSpacing) }; std::vector intPnts; - myHypo->SetGridSpacing(grdSpace, intPnts, 0 ); // Spacing in dir 0 - myHypo->SetGridSpacing(grdSpace, intPnts, 1 ); // Spacing in dir 1 - myHypo->SetGridSpacing(grdSpace, intPnts, 2 ); // Spacing in dir 2 - myHypo->SetSizeThreshold(theSizeThreshold); // set threshold - myHypo->GetCoordinates(xCoords, yCoords, zCoords, shapeBox); - grid.SetCoordinates( xCoords, yCoords, zCoords, axisDirs, shapeBox ); - - for ( size_t i = 0; i < facesItersectors.size(); ++i ) - facesItersectors[i].Intersect(); - - for ( size_t i = 0; i < facesItersectors.size(); ++i ) - facesItersectors[i].StoreIntersections(); - - grid.ComputeNodes( *grid._helper ); - grid.GetEdgesToImplement( edge2faceIDsMap, theShape, faceVec ); + + std::unique_ptr aHypo ( new CartesianHypo() ); + aHypo->SetAxisDirs(axisDirs); + aHypo->SetGridSpacing(grdSpace, intPnts, 0 ); // Spacing in dir 0 + aHypo->SetGridSpacing(grdSpace, intPnts, 1 ); // Spacing in dir 1 + aHypo->SetGridSpacing(grdSpace, intPnts, 2 ); // Spacing in dir 2 + aHypo->SetSizeThreshold(theSizeThreshold); // set threshold + + grid.GridInitAndInterserctWithShape(theShape, edge2faceIDsMap, aHypo.get(), theNumOfThreads, false); } /*! @@ -139,6 +106,10 @@ void GridInitAndIntersectWithShape( Grid& grid, */ bool testNRTM1() { + TopoDS_Shape aShape; + loadBrepShape( "data/HexahedronTest/NRTM1.brep", aShape ); + CPPUNIT_ASSERT_MESSAGE( "Could not load the brep shape!", !aShape.IsNull() ); + const auto numOfCores = std::thread::hardware_concurrency() == 0 ? 1 : std::thread::hardware_concurrency(); std::vector numberOfThreads(numOfCores); std::iota (std::begin(numberOfThreads), std::end(numberOfThreads), 1); @@ -147,20 +118,22 @@ bool testNRTM1() { for (size_t i = 0; i < 10 /*trials*/; i++) { - TopoDS_Shape myShape; - loadBrepShape( "data/HexahedronTest/NRTM1.brep", myShape ); - CPPUNIT_ASSERT_MESSAGE( "Could not load the brep shape!", !myShape.IsNull() ); - std::unique_ptr myMesh( new SMESH_Mesh_Test() ); - myMesh->ShapeToMesh( myShape ); - SMESH_MesherHelper helper( *myMesh ); + std::unique_ptr aMesh( new SMESH_Mesh_Test() ); + aMesh->ShapeToMesh( aShape ); + SMESH_MesherHelper helper( *aMesh ); + Grid grid; grid._helper = &helper; - grid._toAddEdges = false; grid._toCreateFaces = false; grid._toConsiderInternalFaces = false; grid._toUseThresholdForInternalFaces = false; grid._toUseQuanta = false; + grid._toAddEdges = false; + grid._toCreateFaces = false; + grid._toConsiderInternalFaces = false; + grid._toUseThresholdForInternalFaces = false; + grid._toUseQuanta = false; grid._sizeThreshold = 4.0; - grid.InitGeometry( myShape ); - std::map< TGeomID, vector< TGeomID > > edge2faceIDsMap; - GridInitAndIntersectWithShape( grid, 1.0, 4.0, myShape, edge2faceIDsMap, nThreads ); + TEdge2faceIDsMap edge2faceIDsMap; + GridInitAndIntersectWithShape( grid, 1.0, 4.0, aShape, edge2faceIDsMap, nThreads ); + Hexahedron hex( &grid ); int nbAdded = hex.MakeElements( helper, edge2faceIDsMap, nThreads ); CPPUNIT_ASSERT_MESSAGE( "Number of computed elements does not match", nbAdded == 1024 ); @@ -174,6 +147,10 @@ bool testNRTM1() */ bool testNRTJ4() { + TopoDS_Shape aShape; + loadBrepShape( "data/HexahedronTest/NRTMJ4.brep", aShape ); + CPPUNIT_ASSERT_MESSAGE( "Could not load the brep shape!", !aShape.IsNull() ); + const auto numOfCores = std::thread::hardware_concurrency() == 0 ? 1 : std::thread::hardware_concurrency()/2; std::vector numberOfThreads(numOfCores); std::iota (std::begin(numberOfThreads), std::end(numberOfThreads), 1); @@ -183,22 +160,22 @@ bool testNRTJ4() { for (size_t i = 0; i < 10 /*trials*/; i++) { - TopoDS_Shape myShape; - loadBrepShape( "data/HexahedronTest/NRTMJ4.brep", myShape ); - CPPUNIT_ASSERT_MESSAGE( "Could not load the brep shape!", !myShape.IsNull() ); - std::unique_ptr myMesh( new SMESH_Mesh_Test() ); - myMesh->ShapeToMesh( myShape ); - SMESH_MesherHelper helper( *myMesh ); + std::unique_ptr aMesh( new SMESH_Mesh_Test() ); + aMesh->ShapeToMesh( aShape ); + SMESH_MesherHelper helper( *aMesh ); + Grid grid; grid._helper = &helper; - grid._toAddEdges = false; grid._toConsiderInternalFaces = false; grid._toUseThresholdForInternalFaces = false; grid._toUseQuanta = false; + grid._toAddEdges = false; + grid._toConsiderInternalFaces = false; + grid._toUseThresholdForInternalFaces = false; + grid._toUseQuanta = false; double testThreshold = 1.000001; grid._toCreateFaces = true; grid._sizeThreshold = testThreshold; - grid.InitGeometry( myShape ); - std::map< TGeomID, vector< TGeomID > > edge2faceIDsMap; - GridInitAndIntersectWithShape( grid, 2.0, testThreshold, myShape, edge2faceIDsMap, nThreads ); + TEdge2faceIDsMap edge2faceIDsMap; + GridInitAndIntersectWithShape( grid, 2.0, testThreshold, aShape, edge2faceIDsMap, nThreads ); Hexahedron hex( &grid ); int nbAdded = hex.MakeElements( helper, edge2faceIDsMap, nThreads ); CPPUNIT_ASSERT_MESSAGE( "Number of computed elements does not match", nbAdded == 35150 ); @@ -210,6 +187,9 @@ bool testNRTJ4() // Entry point for test int main() { - bool isOK = testNRTM1() && testNRTJ4(); + bool isOK = testNRTM1(); + if (!testNRTJ4()) + isOK = false; + return isOK ? 0 : 1; } diff --git a/src/StdMeshers/StdMeshers_Cartesian_3D.cxx b/src/StdMeshers/StdMeshers_Cartesian_3D.cxx index 98733b703..c20c5a2a0 100644 --- a/src/StdMeshers/StdMeshers_Cartesian_3D.cxx +++ b/src/StdMeshers/StdMeshers_Cartesian_3D.cxx @@ -52,8 +52,9 @@ using namespace std; using namespace SMESH; -using namespace gridtools; +using namespace StdMeshers::Cartesian3D; // using namespace facegridintersector; + namespace { /*! @@ -229,7 +230,7 @@ bool StdMeshers_Cartesian_3D::Compute(SMESH_Mesh & theMesh, } // remove free nodes - if ( SMESHDS_SubMesh * smDS = meshDS->MeshElements( helper.GetSubShapeID() )) + if ( /*SMESHDS_SubMesh * smDS = */meshDS->MeshElements( helper.GetSubShapeID() )) { std::vector< const SMDS_MeshNode* > nodesToRemove; // get intersection nodes @@ -391,4 +392,4 @@ void StdMeshers_Cartesian_3D::setSubmeshesComputed(SMESH_Mesh& theMesh, { for ( TopExp_Explorer soExp( theShape, TopAbs_SOLID ); soExp.More(); soExp.Next() ) _EventListener::setAlwaysComputed( true, theMesh.GetSubMesh( soExp.Current() )); -} \ No newline at end of file +} diff --git a/src/StdMeshers/StdMeshers_Cartesian_3D_Grid.cxx b/src/StdMeshers/StdMeshers_Cartesian_3D_Grid.cxx index 0d02a9f73..273b51c19 100644 --- a/src/StdMeshers/StdMeshers_Cartesian_3D_Grid.cxx +++ b/src/StdMeshers/StdMeshers_Cartesian_3D_Grid.cxx @@ -35,9 +35,12 @@ #include #endif +using namespace std; using namespace SMESH; -using namespace gridtools; +using namespace StdMeshers::Cartesian3D; + std::mutex _bMutex; + //============================================================================= /* * Remove coincident intersection points @@ -929,8 +932,10 @@ void Grid::ComputeNodes(SMESH_MesherHelper& helper) return; } -bool Grid::GridInitAndInterserctWithShape( const TopoDS_Shape& theShape, std::map< TGeomID, vector< TGeomID > >& edge2faceIDsMap, - const StdMeshers_CartesianParameters3D* hyp, const int numOfThreads, bool computeCanceled ) +bool Grid::GridInitAndInterserctWithShape( const TopoDS_Shape& theShape, + TEdge2faceIDsMap& edge2faceIDsMap, + const StdMeshers_CartesianParameters3D* hyp, + const int /*numOfThreads*/, bool computeCanceled ) { InitGeometry( theShape ); @@ -1088,3 +1093,351 @@ void Tools::GetExactBndBox( const vector< TopoDS_Shape >& faceVec, const double* return; } +//============================================================================= +/* + * Intersects TopoDS_Face with all GridLine's + */ +void FaceGridIntersector::Intersect() +{ + FaceLineIntersector intersector; + intersector._surfaceInt = GetCurveFaceIntersector(); + intersector._tol = _grid->_tol; + intersector._transOut = _face.Orientation() == TopAbs_REVERSED ? Trans_IN : Trans_OUT; + intersector._transIn = _face.Orientation() == TopAbs_REVERSED ? Trans_OUT : Trans_IN; + + typedef void (FaceLineIntersector::* PIntFun )(const GridLine& gridLine); + PIntFun interFunction; + + bool isDirect = true; + BRepAdaptor_Surface surf( _face ); + switch ( surf.GetType() ) { + case GeomAbs_Plane: + intersector._plane = surf.Plane(); + interFunction = &FaceLineIntersector::IntersectWithPlane; + isDirect = intersector._plane.Direct(); + break; + case GeomAbs_Cylinder: + intersector._cylinder = surf.Cylinder(); + interFunction = &FaceLineIntersector::IntersectWithCylinder; + isDirect = intersector._cylinder.Direct(); + break; + case GeomAbs_Cone: + intersector._cone = surf.Cone(); + interFunction = &FaceLineIntersector::IntersectWithCone; + //isDirect = intersector._cone.Direct(); + break; + case GeomAbs_Sphere: + intersector._sphere = surf.Sphere(); + interFunction = &FaceLineIntersector::IntersectWithSphere; + isDirect = intersector._sphere.Direct(); + break; + case GeomAbs_Torus: + intersector._torus = surf.Torus(); + interFunction = &FaceLineIntersector::IntersectWithTorus; + //isDirect = intersector._torus.Direct(); + break; + default: + interFunction = &FaceLineIntersector::IntersectWithSurface; + } + if ( !isDirect ) + std::swap( intersector._transOut, intersector._transIn ); + + _intersections.clear(); + for ( int iDir = 0; iDir < 3; ++iDir ) // loop on 3 line directions + { + if ( surf.GetType() == GeomAbs_Plane ) + { + // check if all lines in this direction are parallel to a plane + if ( intersector._plane.Axis().IsNormal( _grid->_lines[iDir][0]._line.Position(), + Precision::Angular())) + continue; + // find out a transition, that is the same for all lines of a direction + gp_Dir plnNorm = intersector._plane.Axis().Direction(); + gp_Dir lineDir = _grid->_lines[iDir][0]._line.Direction(); + intersector._transition = + ( plnNorm * lineDir < 0 ) ? intersector._transIn : intersector._transOut; + } + if ( surf.GetType() == GeomAbs_Cylinder ) + { + // check if all lines in this direction are parallel to a cylinder + if ( intersector._cylinder.Axis().IsParallel( _grid->_lines[iDir][0]._line.Position(), + Precision::Angular())) + continue; + } + + // intersect the grid lines with the face + for ( size_t iL = 0; iL < _grid->_lines[iDir].size(); ++iL ) + { + GridLine& gridLine = _grid->_lines[iDir][iL]; + if ( _bndBox.IsOut( gridLine._line )) continue; + + intersector._intPoints.clear(); + (intersector.*interFunction)( gridLine ); // <- intersection with gridLine + for ( size_t i = 0; i < intersector._intPoints.size(); ++i ) + _intersections.push_back( std::make_pair( &gridLine, intersector._intPoints[i] )); + } + } + + if ( _face.Orientation() == TopAbs_INTERNAL ) + { + for ( size_t i = 0; i < _intersections.size(); ++i ) + if ( _intersections[i].second._transition == Trans_IN || + _intersections[i].second._transition == Trans_OUT ) + { + _intersections[i].second._transition = Trans_INTERNAL; + } + } + return; +} + +#ifdef WITH_TBB +//================================================================================ +/* + * check if its face can be safely intersected in a thread + */ +bool FaceGridIntersector::IsThreadSafe(std::set< const Standard_Transient* >& noSafeTShapes) const +{ + bool isSafe = true; + + // check surface + TopLoc_Location loc; + Handle(Geom_Surface) surf = BRep_Tool::Surface( _face, loc ); + Handle(Geom_RectangularTrimmedSurface) ts = + Handle(Geom_RectangularTrimmedSurface)::DownCast( surf ); + while( !ts.IsNull() ) { + surf = ts->BasisSurface(); + ts = Handle(Geom_RectangularTrimmedSurface)::DownCast(surf); + } + if ( surf->IsKind( STANDARD_TYPE(Geom_BSplineSurface )) || + surf->IsKind( STANDARD_TYPE(Geom_BezierSurface ))) + if ( !noSafeTShapes.insert( _face.TShape().get() ).second ) + isSafe = false; + + double f, l; + TopExp_Explorer exp( _face, TopAbs_EDGE ); + for ( ; exp.More(); exp.Next() ) + { + bool edgeIsSafe = true; + const TopoDS_Edge& e = TopoDS::Edge( exp.Current() ); + // check 3d curve + { + Handle(Geom_Curve) c = BRep_Tool::Curve( e, loc, f, l); + if ( !c.IsNull() ) + { + Handle(Geom_TrimmedCurve) tc = Handle(Geom_TrimmedCurve)::DownCast(c); + while( !tc.IsNull() ) { + c = tc->BasisCurve(); + tc = Handle(Geom_TrimmedCurve)::DownCast(c); + } + if ( c->IsKind( STANDARD_TYPE(Geom_BSplineCurve )) || + c->IsKind( STANDARD_TYPE(Geom_BezierCurve ))) + edgeIsSafe = false; + } + } + // check 2d curve + if ( edgeIsSafe ) + { + Handle(Geom2d_Curve) c2 = BRep_Tool::CurveOnSurface( e, surf, loc, f, l); + if ( !c2.IsNull() ) + { + Handle(Geom2d_TrimmedCurve) tc = Handle(Geom2d_TrimmedCurve)::DownCast(c2); + while( !tc.IsNull() ) { + c2 = tc->BasisCurve(); + tc = Handle(Geom2d_TrimmedCurve)::DownCast(c2); + } + if ( c2->IsKind( STANDARD_TYPE(Geom2d_BSplineCurve )) || + c2->IsKind( STANDARD_TYPE(Geom2d_BezierCurve ))) + edgeIsSafe = false; + } + } + if ( !edgeIsSafe && !noSafeTShapes.insert( e.TShape().get() ).second ) + isSafe = false; + } + return isSafe; +} +#endif + +//================================================================================ +/* + * Store an intersection if it is IN or ON the face + */ +void FaceLineIntersector::addIntPoint(const bool toClassify) +{ + if ( !toClassify || UVIsOnFace() ) + { + F_IntersectPoint p; + p._paramOnLine = _w; + p._u = _u; + p._v = _v; + p._transition = _transition; + _intPoints.push_back( p ); + } +} + +//================================================================================ +/* + * Intersect a line with a plane + */ +void FaceLineIntersector::IntersectWithPlane(const GridLine& gridLine) +{ + IntAna_IntConicQuad linPlane( gridLine._line, _plane, Precision::Angular()); + _w = linPlane.ParamOnConic(1); + if ( isParamOnLineOK( gridLine._length )) + { + ElSLib::Parameters(_plane, linPlane.Point(1) ,_u,_v); + addIntPoint(); + } +} + +//================================================================================ +/* + * Intersect a line with a cylinder + */ +void FaceLineIntersector::IntersectWithCylinder(const GridLine& gridLine) +{ + IntAna_IntConicQuad linCylinder( gridLine._line, _cylinder ); + if ( linCylinder.IsDone() && linCylinder.NbPoints() > 0 ) + { + _w = linCylinder.ParamOnConic(1); + if ( linCylinder.NbPoints() == 1 ) + _transition = Trans_TANGENT; + else + _transition = _w < linCylinder.ParamOnConic(2) ? _transIn : _transOut; + if ( isParamOnLineOK( gridLine._length )) + { + ElSLib::Parameters(_cylinder, linCylinder.Point(1) ,_u,_v); + addIntPoint(); + } + if ( linCylinder.NbPoints() > 1 ) + { + _w = linCylinder.ParamOnConic(2); + if ( isParamOnLineOK( gridLine._length )) + { + ElSLib::Parameters(_cylinder, linCylinder.Point(2) ,_u,_v); + _transition = ( _transition == Trans_OUT ) ? Trans_IN : Trans_OUT; + addIntPoint(); + } + } + } +} + +//================================================================================ +/* + * Intersect a line with a cone + */ +void FaceLineIntersector::IntersectWithCone (const GridLine& gridLine) +{ + IntAna_IntConicQuad linCone(gridLine._line,_cone); + if ( !linCone.IsDone() ) return; + gp_Pnt P; + gp_Vec du, dv, norm; + for ( int i = 1; i <= linCone.NbPoints(); ++i ) + { + _w = linCone.ParamOnConic( i ); + if ( !isParamOnLineOK( gridLine._length )) continue; + ElSLib::Parameters(_cone, linCone.Point(i) ,_u,_v); + if ( UVIsOnFace() ) + { + ElSLib::D1( _u, _v, _cone, P, du, dv ); + norm = du ^ dv; + double normSize2 = norm.SquareMagnitude(); + if ( normSize2 > Precision::Angular() * Precision::Angular() ) + { + double cos = norm.XYZ() * gridLine._line.Direction().XYZ(); + cos /= sqrt( normSize2 ); + if ( cos < -Precision::Angular() ) + _transition = _transIn; + else if ( cos > Precision::Angular() ) + _transition = _transOut; + else + _transition = Trans_TANGENT; + } + else + { + _transition = Trans_APEX; + } + addIntPoint( /*toClassify=*/false); + } + } +} + +//================================================================================ +/* + * Intersect a line with a sphere + */ +void FaceLineIntersector::IntersectWithSphere (const GridLine& gridLine) +{ + IntAna_IntConicQuad linSphere(gridLine._line,_sphere); + if ( linSphere.IsDone() && linSphere.NbPoints() > 0 ) + { + _w = linSphere.ParamOnConic(1); + if ( linSphere.NbPoints() == 1 ) + _transition = Trans_TANGENT; + else + _transition = _w < linSphere.ParamOnConic(2) ? _transIn : _transOut; + if ( isParamOnLineOK( gridLine._length )) + { + ElSLib::Parameters(_sphere, linSphere.Point(1) ,_u,_v); + addIntPoint(); + } + if ( linSphere.NbPoints() > 1 ) + { + _w = linSphere.ParamOnConic(2); + if ( isParamOnLineOK( gridLine._length )) + { + ElSLib::Parameters(_sphere, linSphere.Point(2) ,_u,_v); + _transition = ( _transition == Trans_OUT ) ? Trans_IN : Trans_OUT; + addIntPoint(); + } + } + } +} + +//================================================================================ +/* + * Intersect a line with a torus + */ +void FaceLineIntersector::IntersectWithTorus (const GridLine& gridLine) +{ + IntAna_IntLinTorus linTorus(gridLine._line,_torus); + if ( !linTorus.IsDone()) return; + gp_Pnt P; + gp_Vec du, dv, norm; + for ( int i = 1; i <= linTorus.NbPoints(); ++i ) + { + _w = linTorus.ParamOnLine( i ); + if ( !isParamOnLineOK( gridLine._length )) continue; + linTorus.ParamOnTorus( i, _u,_v ); + if ( UVIsOnFace() ) + { + ElSLib::D1( _u, _v, _torus, P, du, dv ); + norm = du ^ dv; + double normSize = norm.Magnitude(); + double cos = norm.XYZ() * gridLine._line.Direction().XYZ(); + cos /= normSize; + if ( cos < -Precision::Angular() ) + _transition = _transIn; + else if ( cos > Precision::Angular() ) + _transition = _transOut; + else + _transition = Trans_TANGENT; + addIntPoint( /*toClassify=*/false); + } + } +} + +//================================================================================ +/* + * Intersect a line with a non-analytical surface + */ +void FaceLineIntersector::IntersectWithSurface (const GridLine& gridLine) +{ + _surfaceInt->Perform( gridLine._line, 0.0, gridLine._length ); + if ( !_surfaceInt->IsDone() ) return; + for ( int i = 1; i <= _surfaceInt->NbPnt(); ++i ) + { + _transition = Transition( _surfaceInt->Transition( i ) ); + _w = _surfaceInt->WParameter( i ); + addIntPoint(/*toClassify=*/false); + } +} diff --git a/src/StdMeshers/StdMeshers_Cartesian_3D_Grid.hxx b/src/StdMeshers/StdMeshers_Cartesian_3D_Grid.hxx index 3f8cc6800..3cf439aa3 100644 --- a/src/StdMeshers/StdMeshers_Cartesian_3D_Grid.hxx +++ b/src/StdMeshers/StdMeshers_Cartesian_3D_Grid.hxx @@ -32,12 +32,12 @@ // STD #include +#include #include #include #include #include - // SMESH #include "SMESH_StdMeshers.hxx" #include "StdMeshers_FaceSide.hxx" @@ -116,14 +116,16 @@ // All utility structs used in Grid and hexahedron class will be included here // Ideally each one of this should define their own testable class -using namespace std; -using namespace SMESH; -namespace gridtools +namespace StdMeshers +{ +namespace Cartesian3D { typedef int TGeomID; // IDs of sub-shapes typedef TopTools_ShapeMapHasher TShapeHasher; // non-oriented shape hasher typedef std::array< int, 3 > TIJK; + typedef std::map< TGeomID, std::vector< TGeomID > > TEdge2faceIDsMap; + const TGeomID theUndefID = 1e+9; //============================================================================= @@ -136,6 +138,7 @@ namespace gridtools Trans_APEX, Trans_INTERNAL // for INTERNAL FACE }; + // -------------------------------------------------------------------------- /*! * \brief Sub-entities of a FACE neighboring its concave VERTEX. @@ -155,6 +158,7 @@ namespace gridtools void SetVertex( TGeomID v ) { ( _v1 ? _v2 : _v1 ) = v; } }; typedef NCollection_DataMap< TGeomID, ConcaveFace > TConcaveVertex2Face; + // -------------------------------------------------------------------------- /*! * \brief Container of IDs of SOLID sub-shapes @@ -167,7 +171,7 @@ namespace gridtools public: virtual ~Solid() {} virtual bool Contains( TGeomID /*subID*/ ) const { return true; } - virtual bool ContainsAny( const vector< TGeomID>& /*subIDs*/ ) const { return true; } + virtual bool ContainsAny( const std::vector< TGeomID>& /*subIDs*/ ) const { return true; } virtual TopAbs_Orientation Orientation( const TopoDS_Shape& s ) const { return s.Orientation(); } virtual bool IsOutsideOriented( TGeomID /*faceID*/ ) const { return true; } void SetID( TGeomID id ) { _id = id; } @@ -179,6 +183,7 @@ namespace gridtools bool HasConcaveVertex() const { return !_concaveVertex.IsEmpty(); } const ConcaveFace* GetConcave( TGeomID V ) const { return _concaveVertex.Seek( V ); } }; + // -------------------------------------------------------------------------- class OneOfSolids : public Solid { @@ -190,7 +195,7 @@ namespace gridtools TopAbs_ShapeEnum subType, const SMESHDS_Mesh* mesh ); virtual bool Contains( TGeomID i ) const { return i == ID() || _subIDs.Contains( i ); } - virtual bool ContainsAny( const vector< TGeomID>& subIDs ) const + virtual bool ContainsAny( const std::vector< TGeomID>& subIDs ) const { for ( size_t i = 0; i < subIDs.size(); ++i ) if ( Contains( subIDs[ i ])) return true; return false; @@ -205,6 +210,7 @@ namespace gridtools return faceID == 0 || _outFaceIDs.Contains( faceID ); } }; + // -------------------------------------------------------------------------- /*! * \brief Hold a vector of TGeomID and clear it at destruction @@ -246,6 +252,7 @@ namespace gridtools return common; } }; + // -------------------------------------------------------------------------- /*! * \brief Geom data @@ -253,20 +260,20 @@ namespace gridtools struct Geometry { TopoDS_Shape _mainShape; - vector< vector< TGeomID > > _solidIDsByShapeID;// V/E/F ID -> SOLID IDs - Solid _soleSolid; - map< TGeomID, OneOfSolids > _solidByID; + std::vector< std::vector< TGeomID > > _solidIDsByShapeID;// V/E/F ID -> SOLID IDs + Solid _soleSolid; + std::map< TGeomID, OneOfSolids > _solidByID; TColStd_MapOfInteger _boundaryFaces; // FACEs on boundary of mesh->ShapeToMesh() TColStd_MapOfInteger _strangeEdges; // EDGEs shared by strange FACEs TGeomID _extIntFaceID; // pseudo FACE - extension of INTERNAL FACE TopTools_DataMapOfShapeInteger _shape2NbNodes; // nb of pre-existing nodes on shapes - Controls::ElementsOnShape _edgeClassifier; - Controls::ElementsOnShape _vertexClassifier; + SMESH::Controls::ElementsOnShape _edgeClassifier; + SMESH::Controls::ElementsOnShape _vertexClassifier; bool IsOneSolid() const { return _solidByID.size() < 2; } - GeomIDVecHelder GetSolidIDsByShapeID( const vector< TGeomID >& shapeIDs ) const; + GeomIDVecHelder GetSolidIDsByShapeID( const std::vector< TGeomID >& shapeIDs ) const; }; // -------------------------------------------------------------------------- @@ -279,16 +286,17 @@ namespace gridtools // See Add method modify _node and _faceIDs class members dinamicaly during execution // of Hexahedron.compute() method. // std::mutex _mutex; - mutable const SMDS_MeshNode* _node; - mutable vector< TGeomID > _faceIDs; + mutable const SMDS_MeshNode* _node; + mutable std::vector< TGeomID > _faceIDs; B_IntersectPoint(): _node(NULL) {} - bool Add( const vector< TGeomID >& fIDs, const SMDS_MeshNode* n=NULL ) const; + bool Add( const std::vector< TGeomID >& fIDs, const SMDS_MeshNode* n=NULL ) const; TGeomID HasCommonFace( const B_IntersectPoint * other, TGeomID avoidFace=-1 ) const; size_t GetCommonFaces( const B_IntersectPoint * other, TGeomID * commonFaces ) const; bool IsOnFace( TGeomID faceID ) const; virtual ~B_IntersectPoint() {} }; + // -------------------------------------------------------------------------- /*! * \brief Data of intersection between a GridLine and a TopoDS_Face @@ -304,6 +312,7 @@ namespace gridtools return _paramOnLine < o._paramOnLine; } }; + // -------------------------------------------------------------------------- /*! * \brief Data of intersection between GridPlanes and a TopoDS_EDGE @@ -323,10 +332,10 @@ namespace gridtools { gp_Lin _line; double _length; // line length - multiset< F_IntersectPoint > _intPoints; + std::multiset< F_IntersectPoint > _intPoints; void RemoveExcessIntPoints( const double tol ); - TGeomID GetSolidIDBefore( multiset< F_IntersectPoint >::iterator ip, + TGeomID GetSolidIDBefore( std::multiset< F_IntersectPoint >::iterator ip, const TGeomID prevID, const Geometry& geom); }; @@ -336,9 +345,9 @@ namespace gridtools */ struct GridPlanes { - gp_XYZ _zNorm; - vector< gp_XYZ > _origins; // origin points of all planes in one direction - vector< double > _zProjs; // projections of origins to _zNorm + gp_XYZ _zNorm; + std::vector< gp_XYZ > _origins; // origin points of all planes in one direction + std::vector< double > _zProjs; // projections of origins to _zNorm }; // -------------------------------------------------------------------------- /*! @@ -349,11 +358,11 @@ namespace gridtools size_t _size [3]; size_t _curInd[3]; size_t _iVar1, _iVar2, _iConst; - string _name1, _name2, _nameConst; + std::string _name1, _name2, _nameConst; LineIndexer() {} LineIndexer( size_t sz1, size_t sz2, size_t sz3, size_t iv1, size_t iv2, size_t iConst, - const string& nv1, const string& nv2, const string& nConst ) + const std::string& nv1, const std::string& nv2, const std::string& nConst ) { _size[0] = sz1; _size[1] = sz2; _size[2] = sz3; _curInd[0] = _curInd[1] = _curInd[2] = 0; @@ -398,31 +407,28 @@ namespace gridtools * \brief computes exact bounding box with axes parallel to given ones */ //================================================================================ - static void GetExactBndBox( const vector< TopoDS_Shape >& faceVec, const double* axesDirs, Bnd_Box& shapeBox ); + static void GetExactBndBox( const std::vector< TopoDS_Shape >& faceVec, const double* axesDirs, Bnd_Box& shapeBox ); }; -} // end namespace gridtools - -using namespace gridtools; -class STDMESHERS_EXPORT Grid -{ -public: - vector< double > _coords[3]; // coordinates of grid nodes - gp_XYZ _axes [3]; // axis directions - vector< GridLine > _lines [3]; // in 3 directions - double _tol, _minCellSize; - gp_XYZ _origin; - gp_Mat _invB; // inverted basis of _axes + class STDMESHERS_EXPORT Grid + { + public: + std::vector< double > _coords[3]; // coordinates of grid nodes + gp_XYZ _axes [3]; // axis directions + std::vector< GridLine > _lines [3]; // in 3 directions + double _tol, _minCellSize; + gp_XYZ _origin; + gp_Mat _invB; // inverted basis of _axes // index shift within _nodes of nodes of a cell from the 1st node - int _nodeShift[8]; + int _nodeShift[8]; - vector< const SMDS_MeshNode* > _nodes; // mesh nodes at grid nodes - vector< const SMDS_MeshNode* > _allBorderNodes; // mesh nodes between the bounding box and the geometry boundary + std::vector< const SMDS_MeshNode* > _nodes; // mesh nodes at grid nodes + std::vector< const SMDS_MeshNode* > _allBorderNodes; // mesh nodes between the bounding box and the geometry boundary - vector< const F_IntersectPoint* > _gridIntP; // grid node intersection with geometry - ObjectPool< E_IntersectPoint > _edgeIntPool; // intersections with EDGEs - ObjectPool< F_IntersectPoint > _extIntPool; // intersections with extended INTERNAL FACEs + std::vector< const F_IntersectPoint* > _gridIntP; // grid node intersection with geometry + ObjectPool< E_IntersectPoint > _edgeIntPool; // intersections with EDGEs + ObjectPool< F_IntersectPoint > _extIntPool; // intersections with extended INTERNAL FACEs //list< E_IntersectPoint > _edgeIntP; // intersections with EDGEs Geometry _geometry; @@ -469,10 +475,10 @@ public: void InitGeometry( const TopoDS_Shape& theShape ); void InitClassifier( const TopoDS_Shape& mainShape, TopAbs_ShapeEnum shapeType, - Controls::ElementsOnShape& classifier ); - void GetEdgesToImplement( map< TGeomID, vector< TGeomID > > & edge2faceMap, + SMESH::Controls::ElementsOnShape& classifier ); + void GetEdgesToImplement( std::map< TGeomID, std::vector< TGeomID > > & edge2faceMap, const TopoDS_Shape& shape, - const vector< TopoDS_Shape >& faces ); + const std::vector< TopoDS_Shape >& faces ); void SetSolidFather( const TopoDS_Shape& s, const TopoDS_Shape& theShapeToMesh ); bool IsShared( TGeomID faceID ) const; bool IsAnyShared( const std::vector< TGeomID >& faceIDs ) const; @@ -486,7 +492,7 @@ public: TGeomID PseudoIntExtFaceID() const { return _geometry._extIntFaceID; } Solid* GetSolid( TGeomID solidID = 0 ); Solid* GetOneOfSolids( TGeomID solidID ); - const vector< TGeomID > & GetSolidIDs( TGeomID subShapeID ) const; + const std::vector< TGeomID > & GetSolidIDs( TGeomID subShapeID ) const; bool IsCorrectTransition( TGeomID faceID, const Solid* solid ); bool IsBoundaryFace( TGeomID face ) const { return _geometry._boundaryFaces.Contains( face ); } void SetOnShape( const SMDS_MeshNode* n, const F_IntersectPoint& ip, @@ -495,19 +501,19 @@ public: bool IsToCheckNodePos() const { return !_toAddEdges && _toCreateFaces; } bool IsToRemoveExcessEntities() const { return !_toAddEdges; } - void SetCoordinates(const vector& xCoords, - const vector& yCoords, - const vector& zCoords, + void SetCoordinates(const std::vector& xCoords, + const std::vector& yCoords, + const std::vector& zCoords, const double* axesDirs, const Bnd_Box& bndBox ); void ComputeUVW(const gp_XYZ& p, double uvw[3]); void ComputeNodes(SMESH_MesherHelper& helper); - bool GridInitAndInterserctWithShape( const TopoDS_Shape& theShape, std::map< TGeomID, vector< TGeomID > >& edge2faceIDsMap, - const StdMeshers_CartesianParameters3D* hyp, const int numOfThreads, bool computeCanceled ); -}; - -namespace -{ + bool GridInitAndInterserctWithShape( const TopoDS_Shape& theShape, + std::map< TGeomID, std::vector< TGeomID > >& edge2faceIDsMap, + const StdMeshers_CartesianParameters3D* hyp, + const int numOfThreads, + bool computeCanceled ); + }; // Implement parallel computation of Hexa with c++ thread implementation template @@ -529,6 +535,7 @@ namespace std::for_each(it, last, f); // last steps while we wait for other threads std::for_each(threads.begin(), threads.end(), [](std::thread& x){x.join();}); } + // -------------------------------------------------------------------------- /*! * \brief Intersector of TopoDS_Face with all GridLine's @@ -540,7 +547,7 @@ namespace Grid* _grid; Bnd_Box _bndBox; IntCurvesFace_Intersector* _surfaceInt; - vector< std::pair< GridLine*, F_IntersectPoint > > _intersections; + std::vector< std::pair< GridLine*, F_IntersectPoint > > _intersections; FaceGridIntersector(): _grid(0), _surfaceInt(0) {} void Intersect(); @@ -549,7 +556,7 @@ namespace { for ( size_t i = 0; i < _intersections.size(); ++i ) { - multiset< F_IntersectPoint >::iterator ip = + std::multiset< F_IntersectPoint >::iterator ip = _intersections[i].first->_intPoints.insert( _intersections[i].second ); ip->_faceIDs.reserve( 1 ); ip->_faceIDs.push_back( _faceID ); @@ -572,7 +579,7 @@ namespace return _surfaceInt; } #ifdef WITH_TBB - bool IsThreadSafe(set< const Standard_Transient* >& noSafeTShapes) const; + bool IsThreadSafe(std::set< const Standard_Transient* >& noSafeTShapes) const; #endif }; @@ -583,8 +590,8 @@ namespace */ struct ParallelIntersector { - vector< FaceGridIntersector >& _faceVec; - ParallelIntersector( vector< FaceGridIntersector >& faceVec): _faceVec(faceVec){} + std::vector< FaceGridIntersector >& _faceVec; + ParallelIntersector( std::vector< FaceGridIntersector >& faceVec): _faceVec(faceVec){} void operator() ( const tbb::blocked_range& r ) const { for ( size_t i = r.begin(); i != r.end(); ++i ) @@ -617,7 +624,7 @@ namespace gp_Torus _torus; IntCurvesFace_Intersector* _surfaceInt; - vector< F_IntersectPoint > _intPoints; + std::vector< F_IntersectPoint > _intPoints; void IntersectWithPlane (const GridLine& gridLine); void IntersectWithCylinder(const GridLine& gridLine); @@ -626,7 +633,14 @@ namespace void IntersectWithTorus (const GridLine& gridLine); void IntersectWithSurface (const GridLine& gridLine); - bool UVIsOnFace() const; + /* + * Return true if (_u,_v) is on the face + */ + bool UVIsOnFace() const + { + TopAbs_State state = _surfaceInt->ClassifyUVPoint(gp_Pnt2d( _u,_v )); + return ( state == TopAbs_IN || state == TopAbs_ON ); + } void addIntPoint(const bool toClassify=true); bool isParamOnLineOK( const double linLength ) { @@ -636,356 +650,7 @@ namespace ~FaceLineIntersector() { if (_surfaceInt ) delete _surfaceInt; _surfaceInt = 0; } }; - //============================================================================= - /* - * Intersects TopoDS_Face with all GridLine's - */ - void FaceGridIntersector::Intersect() - { - FaceLineIntersector intersector; - intersector._surfaceInt = GetCurveFaceIntersector(); - intersector._tol = _grid->_tol; - intersector._transOut = _face.Orientation() == TopAbs_REVERSED ? Trans_IN : Trans_OUT; - intersector._transIn = _face.Orientation() == TopAbs_REVERSED ? Trans_OUT : Trans_IN; - - typedef void (FaceLineIntersector::* PIntFun )(const GridLine& gridLine); - PIntFun interFunction; - - bool isDirect = true; - BRepAdaptor_Surface surf( _face ); - switch ( surf.GetType() ) { - case GeomAbs_Plane: - intersector._plane = surf.Plane(); - interFunction = &FaceLineIntersector::IntersectWithPlane; - isDirect = intersector._plane.Direct(); - break; - case GeomAbs_Cylinder: - intersector._cylinder = surf.Cylinder(); - interFunction = &FaceLineIntersector::IntersectWithCylinder; - isDirect = intersector._cylinder.Direct(); - break; - case GeomAbs_Cone: - intersector._cone = surf.Cone(); - interFunction = &FaceLineIntersector::IntersectWithCone; - //isDirect = intersector._cone.Direct(); - break; - case GeomAbs_Sphere: - intersector._sphere = surf.Sphere(); - interFunction = &FaceLineIntersector::IntersectWithSphere; - isDirect = intersector._sphere.Direct(); - break; - case GeomAbs_Torus: - intersector._torus = surf.Torus(); - interFunction = &FaceLineIntersector::IntersectWithTorus; - //isDirect = intersector._torus.Direct(); - break; - default: - interFunction = &FaceLineIntersector::IntersectWithSurface; - } - if ( !isDirect ) - std::swap( intersector._transOut, intersector._transIn ); - - _intersections.clear(); - for ( int iDir = 0; iDir < 3; ++iDir ) // loop on 3 line directions - { - if ( surf.GetType() == GeomAbs_Plane ) - { - // check if all lines in this direction are parallel to a plane - if ( intersector._plane.Axis().IsNormal( _grid->_lines[iDir][0]._line.Position(), - Precision::Angular())) - continue; - // find out a transition, that is the same for all lines of a direction - gp_Dir plnNorm = intersector._plane.Axis().Direction(); - gp_Dir lineDir = _grid->_lines[iDir][0]._line.Direction(); - intersector._transition = - ( plnNorm * lineDir < 0 ) ? intersector._transIn : intersector._transOut; - } - if ( surf.GetType() == GeomAbs_Cylinder ) - { - // check if all lines in this direction are parallel to a cylinder - if ( intersector._cylinder.Axis().IsParallel( _grid->_lines[iDir][0]._line.Position(), - Precision::Angular())) - continue; - } - - // intersect the grid lines with the face - for ( size_t iL = 0; iL < _grid->_lines[iDir].size(); ++iL ) - { - GridLine& gridLine = _grid->_lines[iDir][iL]; - if ( _bndBox.IsOut( gridLine._line )) continue; - - intersector._intPoints.clear(); - (intersector.*interFunction)( gridLine ); // <- intersection with gridLine - for ( size_t i = 0; i < intersector._intPoints.size(); ++i ) - _intersections.push_back( make_pair( &gridLine, intersector._intPoints[i] )); - } - } - - if ( _face.Orientation() == TopAbs_INTERNAL ) - { - for ( size_t i = 0; i < _intersections.size(); ++i ) - if ( _intersections[i].second._transition == Trans_IN || - _intersections[i].second._transition == Trans_OUT ) - { - _intersections[i].second._transition = Trans_INTERNAL; - } - } - return; - } - //================================================================================ - /* - * Return true if (_u,_v) is on the face - */ - bool FaceLineIntersector::UVIsOnFace() const - { - TopAbs_State state = _surfaceInt->ClassifyUVPoint(gp_Pnt2d( _u,_v )); - return ( state == TopAbs_IN || state == TopAbs_ON ); - } - //================================================================================ - /* - * Store an intersection if it is IN or ON the face - */ - void FaceLineIntersector::addIntPoint(const bool toClassify) - { - if ( !toClassify || UVIsOnFace() ) - { - F_IntersectPoint p; - p._paramOnLine = _w; - p._u = _u; - p._v = _v; - p._transition = _transition; - _intPoints.push_back( p ); - } - } - //================================================================================ - /* - * Intersect a line with a plane - */ - void FaceLineIntersector::IntersectWithPlane(const GridLine& gridLine) - { - IntAna_IntConicQuad linPlane( gridLine._line, _plane, Precision::Angular()); - _w = linPlane.ParamOnConic(1); - if ( isParamOnLineOK( gridLine._length )) - { - ElSLib::Parameters(_plane, linPlane.Point(1) ,_u,_v); - addIntPoint(); - } - } - //================================================================================ - /* - * Intersect a line with a cylinder - */ - void FaceLineIntersector::IntersectWithCylinder(const GridLine& gridLine) - { - IntAna_IntConicQuad linCylinder( gridLine._line, _cylinder ); - if ( linCylinder.IsDone() && linCylinder.NbPoints() > 0 ) - { - _w = linCylinder.ParamOnConic(1); - if ( linCylinder.NbPoints() == 1 ) - _transition = Trans_TANGENT; - else - _transition = _w < linCylinder.ParamOnConic(2) ? _transIn : _transOut; - if ( isParamOnLineOK( gridLine._length )) - { - ElSLib::Parameters(_cylinder, linCylinder.Point(1) ,_u,_v); - addIntPoint(); - } - if ( linCylinder.NbPoints() > 1 ) - { - _w = linCylinder.ParamOnConic(2); - if ( isParamOnLineOK( gridLine._length )) - { - ElSLib::Parameters(_cylinder, linCylinder.Point(2) ,_u,_v); - _transition = ( _transition == Trans_OUT ) ? Trans_IN : Trans_OUT; - addIntPoint(); - } - } - } - } - //================================================================================ - /* - * Intersect a line with a cone - */ - void FaceLineIntersector::IntersectWithCone (const GridLine& gridLine) - { - IntAna_IntConicQuad linCone(gridLine._line,_cone); - if ( !linCone.IsDone() ) return; - gp_Pnt P; - gp_Vec du, dv, norm; - for ( int i = 1; i <= linCone.NbPoints(); ++i ) - { - _w = linCone.ParamOnConic( i ); - if ( !isParamOnLineOK( gridLine._length )) continue; - ElSLib::Parameters(_cone, linCone.Point(i) ,_u,_v); - if ( UVIsOnFace() ) - { - ElSLib::D1( _u, _v, _cone, P, du, dv ); - norm = du ^ dv; - double normSize2 = norm.SquareMagnitude(); - if ( normSize2 > Precision::Angular() * Precision::Angular() ) - { - double cos = norm.XYZ() * gridLine._line.Direction().XYZ(); - cos /= sqrt( normSize2 ); - if ( cos < -Precision::Angular() ) - _transition = _transIn; - else if ( cos > Precision::Angular() ) - _transition = _transOut; - else - _transition = Trans_TANGENT; - } - else - { - _transition = Trans_APEX; - } - addIntPoint( /*toClassify=*/false); - } - } - } - //================================================================================ - /* - * Intersect a line with a sphere - */ - void FaceLineIntersector::IntersectWithSphere (const GridLine& gridLine) - { - IntAna_IntConicQuad linSphere(gridLine._line,_sphere); - if ( linSphere.IsDone() && linSphere.NbPoints() > 0 ) - { - _w = linSphere.ParamOnConic(1); - if ( linSphere.NbPoints() == 1 ) - _transition = Trans_TANGENT; - else - _transition = _w < linSphere.ParamOnConic(2) ? _transIn : _transOut; - if ( isParamOnLineOK( gridLine._length )) - { - ElSLib::Parameters(_sphere, linSphere.Point(1) ,_u,_v); - addIntPoint(); - } - if ( linSphere.NbPoints() > 1 ) - { - _w = linSphere.ParamOnConic(2); - if ( isParamOnLineOK( gridLine._length )) - { - ElSLib::Parameters(_sphere, linSphere.Point(2) ,_u,_v); - _transition = ( _transition == Trans_OUT ) ? Trans_IN : Trans_OUT; - addIntPoint(); - } - } - } - } - //================================================================================ - /* - * Intersect a line with a torus - */ - void FaceLineIntersector::IntersectWithTorus (const GridLine& gridLine) - { - IntAna_IntLinTorus linTorus(gridLine._line,_torus); - if ( !linTorus.IsDone()) return; - gp_Pnt P; - gp_Vec du, dv, norm; - for ( int i = 1; i <= linTorus.NbPoints(); ++i ) - { - _w = linTorus.ParamOnLine( i ); - if ( !isParamOnLineOK( gridLine._length )) continue; - linTorus.ParamOnTorus( i, _u,_v ); - if ( UVIsOnFace() ) - { - ElSLib::D1( _u, _v, _torus, P, du, dv ); - norm = du ^ dv; - double normSize = norm.Magnitude(); - double cos = norm.XYZ() * gridLine._line.Direction().XYZ(); - cos /= normSize; - if ( cos < -Precision::Angular() ) - _transition = _transIn; - else if ( cos > Precision::Angular() ) - _transition = _transOut; - else - _transition = Trans_TANGENT; - addIntPoint( /*toClassify=*/false); - } - } - } - //================================================================================ - /* - * Intersect a line with a non-analytical surface - */ - void FaceLineIntersector::IntersectWithSurface (const GridLine& gridLine) - { - _surfaceInt->Perform( gridLine._line, 0.0, gridLine._length ); - if ( !_surfaceInt->IsDone() ) return; - for ( int i = 1; i <= _surfaceInt->NbPnt(); ++i ) - { - _transition = Transition( _surfaceInt->Transition( i ) ); - _w = _surfaceInt->WParameter( i ); - addIntPoint(/*toClassify=*/false); - } - } - -#ifdef WITH_TBB - //================================================================================ - /* - * check if its face can be safely intersected in a thread - */ - bool FaceGridIntersector::IsThreadSafe(set< const Standard_Transient* >& noSafeTShapes) const - { - bool isSafe = true; - - // check surface - TopLoc_Location loc; - Handle(Geom_Surface) surf = BRep_Tool::Surface( _face, loc ); - Handle(Geom_RectangularTrimmedSurface) ts = - Handle(Geom_RectangularTrimmedSurface)::DownCast( surf ); - while( !ts.IsNull() ) { - surf = ts->BasisSurface(); - ts = Handle(Geom_RectangularTrimmedSurface)::DownCast(surf); - } - if ( surf->IsKind( STANDARD_TYPE(Geom_BSplineSurface )) || - surf->IsKind( STANDARD_TYPE(Geom_BezierSurface ))) - if ( !noSafeTShapes.insert( _face.TShape().get() ).second ) - isSafe = false; - - double f, l; - TopExp_Explorer exp( _face, TopAbs_EDGE ); - for ( ; exp.More(); exp.Next() ) - { - bool edgeIsSafe = true; - const TopoDS_Edge& e = TopoDS::Edge( exp.Current() ); - // check 3d curve - { - Handle(Geom_Curve) c = BRep_Tool::Curve( e, loc, f, l); - if ( !c.IsNull() ) - { - Handle(Geom_TrimmedCurve) tc = Handle(Geom_TrimmedCurve)::DownCast(c); - while( !tc.IsNull() ) { - c = tc->BasisCurve(); - tc = Handle(Geom_TrimmedCurve)::DownCast(c); - } - if ( c->IsKind( STANDARD_TYPE(Geom_BSplineCurve )) || - c->IsKind( STANDARD_TYPE(Geom_BezierCurve ))) - edgeIsSafe = false; - } - } - // check 2d curve - if ( edgeIsSafe ) - { - Handle(Geom2d_Curve) c2 = BRep_Tool::CurveOnSurface( e, surf, loc, f, l); - if ( !c2.IsNull() ) - { - Handle(Geom2d_TrimmedCurve) tc = Handle(Geom2d_TrimmedCurve)::DownCast(c2); - while( !tc.IsNull() ) { - c2 = tc->BasisCurve(); - tc = Handle(Geom2d_TrimmedCurve)::DownCast(c2); - } - if ( c2->IsKind( STANDARD_TYPE(Geom2d_BSplineCurve )) || - c2->IsKind( STANDARD_TYPE(Geom2d_BezierCurve ))) - edgeIsSafe = false; - } - } - if ( !edgeIsSafe && !noSafeTShapes.insert( e.TShape().get() ).second ) - isSafe = false; - } - return isSafe; - } -#endif -} +} // end namespace Cartesian3D +} // end namespace StdMeshers #endif diff --git a/src/StdMeshers/StdMeshers_Cartesian_3D_Hexahedron.cxx b/src/StdMeshers/StdMeshers_Cartesian_3D_Hexahedron.cxx index f1f8f6965..8a2b10300 100644 --- a/src/StdMeshers/StdMeshers_Cartesian_3D_Hexahedron.cxx +++ b/src/StdMeshers/StdMeshers_Cartesian_3D_Hexahedron.cxx @@ -23,8 +23,9 @@ #include "StdMeshers_Cartesian_3D_Hexahedron.hxx" +using namespace std; using namespace SMESH; -using namespace gridtools; +using namespace StdMeshers::Cartesian3D; std::mutex _eMutex; @@ -856,8 +857,31 @@ void Hexahedron::computeElements( const Solid* solid, int solidIndex ) /*! * \brief Compute mesh volumes resulted from intersection of the Hexahedron */ -void Hexahedron::defineHexahedralFaces( std::vector< _OrientedLink >& splits, std::vector<_Node*>& chainNodes, std::set< TGeomID >& concaveFaces, bool toCheckSideDivision ) +void Hexahedron::defineHexahedralFaces( const Solid* solid, const IsInternalFlag intFlag ) { + std::set< TGeomID > concaveFaces; // to avoid connecting nodes laying on them + + if ( solid->HasConcaveVertex() ) { + for ( const E_IntersectPoint* ip : _eIntPoints ) { + if ( const ConcaveFace* cf = solid->GetConcave( ip->_shapeID )) { + if ( this->hasEdgesAround( cf )) + concaveFaces.insert( cf->_concaveFace ); + } + } + if ( concaveFaces.empty() || concaveFaces.size() * 3 < _eIntPoints.size() ) { + for ( const _Node& hexNode: _hexNodes ) { + if ( hexNode._node && hexNode._intPoint && hexNode._intPoint->_faceIDs.size() >= 3 ) + if ( const ConcaveFace* cf = solid->GetConcave( hexNode._node->GetShapeID() )) + if ( this->hasEdgesAround( cf )) + concaveFaces.insert( cf->_concaveFace ); + } + } + } + + const bool toCheckSideDivision = isImplementEdges() || intFlag & IS_CUT_BY_INTERNAL_FACE; + + std::vector< _OrientedLink > splits; + std::vector< _Node* > chainNodes; for ( int iF = 0; iF < 6; ++iF ) // loop on 6 sides of a hexahedron { _Face& quad = _hexQuads[ iF ] ; @@ -1026,6 +1050,47 @@ void Hexahedron::defineHexahedralFaces( std::vector< _OrientedLink >& splits, st } // } +//================================================================================ +/*! + * \brief get a remaining link to start from, one lying on minimal nb of FACEs + */ +Hexahedron::TFaceOfLink Hexahedron::findStartLink(const vector<_OrientedLink*>& freeLinks, + std::set& usedFaceIDs) +{ + TFaceOfLink faceOfLink( -1, -1 ); + TFaceOfLink facesOfLink[3] = { faceOfLink, faceOfLink, faceOfLink }; + for ( size_t iL = 0; iL < freeLinks.size(); ++iL ) { + if ( freeLinks[ iL ] ) { + std::vector< TGeomID > faces = freeLinks[ iL ]->GetNotUsedFaces( usedFaceIDs ); + if ( faces.size() == 1 ) { + faceOfLink = TFaceOfLink( faces[0], iL ); + if ( !freeLinks[ iL ]->HasEdgeNodes() ) + break; + facesOfLink[0] = faceOfLink; + } + else if ( facesOfLink[0].first < 0 ) { + faceOfLink = TFaceOfLink(( faces.empty() ? -1 : faces[0]), iL ); + facesOfLink[ 1 + faces.empty() ] = faceOfLink; + } + } + } + for ( int i = 0; faceOfLink.first < 0 && i < 3; ++i ) { + faceOfLink = facesOfLink[i]; + } + + if ( faceOfLink.first < 0 ) { // all faces used + for ( size_t iL = 0; iL < freeLinks.size() && faceOfLink.first < 1; ++iL ) { + _OrientedLink* curLinkLocal = freeLinks[ iL ]; + if ( curLinkLocal ) { + faceOfLink.first = curLinkLocal->FirstNode()->IsLinked( curLinkLocal->LastNode()->_intPoint ); + faceOfLink.second = iL; + } + } + usedFaceIDs.clear(); + } + return faceOfLink; +} + //================================================================================ /*! * \brief Compute mesh volumes resulted from intersection of the Hexahedron @@ -1041,36 +1106,10 @@ bool Hexahedron::compute( const Solid* solid, const IsInternalFlag intFlag ) if ( intFlag & IS_CUT_BY_INTERNAL_FACE && !_grid->_toAddEdges ) // Issue #19913 preventVolumesOverlapping(); - std::set< TGeomID > concaveFaces; // to avoid connecting nodes laying on them - - if ( solid->HasConcaveVertex() ) - { - for ( const E_IntersectPoint* ip : _eIntPoints ) - { - if ( const ConcaveFace* cf = solid->GetConcave( ip->_shapeID )) - if ( this->hasEdgesAround( cf )) - concaveFaces.insert( cf->_concaveFace ); - } - if ( concaveFaces.empty() || concaveFaces.size() * 3 < _eIntPoints.size() ) - for ( const _Node& hexNode: _hexNodes ) - { - if ( hexNode._node && hexNode._intPoint && hexNode._intPoint->_faceIDs.size() >= 3 ) - if ( const ConcaveFace* cf = solid->GetConcave( hexNode._node->GetShapeID() )) - if ( this->hasEdgesAround( cf )) - concaveFaces.insert( cf->_concaveFace ); - } - } + const bool hasEdgeIntersections = !_eIntPoints.empty(); // Create polygons from quadrangles - // -------------------------------- - - vector< _OrientedLink > splits; - vector<_Node*> chainNodes; - _Face* coplanarPolyg; - - const bool hasEdgeIntersections = !_eIntPoints.empty(); - const bool toCheckSideDivision = isImplementEdges() || intFlag & IS_CUT_BY_INTERNAL_FACE; - defineHexahedralFaces( splits, chainNodes, concaveFaces, toCheckSideDivision ); + defineHexahedralFaces( solid, intFlag ); // Create polygons closing holes in a polyhedron // ---------------------------------------------- @@ -1080,24 +1119,24 @@ bool Hexahedron::compute( const Solid* solid, const IsInternalFlag intFlag ) _intNodes[ iN ]._usedInFace = 0; // add polygons to their links and mark used nodes - for ( size_t iP = 0; iP < _polygons.size(); ++iP ) - { + for ( size_t iP = 0; iP < _polygons.size(); ++iP ) { _Face& polygon = _polygons[ iP ]; - for ( size_t iL = 0; iL < polygon._links.size(); ++iL ) - { + for ( size_t iL = 0; iL < polygon._links.size(); ++iL ) { polygon._links[ iL ].AddFace( &polygon ); polygon._links[ iL ].FirstNode()->_usedInFace = &polygon; } } + // find free links vector< _OrientedLink* > freeLinks; freeLinks.reserve(20); for ( size_t iP = 0; iP < _polygons.size(); ++iP ) { _Face& polygon = _polygons[ iP ]; - for ( size_t iL = 0; iL < polygon._links.size(); ++iL ) + for ( size_t iL = 0; iL < polygon._links.size(); ++iL ) { if ( polygon._links[ iL ].NbFaces() < 2 ) freeLinks.push_back( & polygon._links[ iL ]); + } } int nbFreeLinks = freeLinks.size(); if ( nbFreeLinks == 1 ) return false; @@ -1123,11 +1162,11 @@ bool Hexahedron::compute( const Solid* solid, const IsInternalFlag intFlag ) } std::set usedFaceIDs; - std::vector< TGeomID > faces; TGeomID curFace = 0; const size_t nbQuadPolygons = _polygons.size(); E_IntersectPoint ipTmp; std::map< TGeomID, std::vector< const B_IntersectPoint* > > tmpAddedFace; // face added to _intPoint + // std::cout << "2\n"; // create polygons by making closed chains of free links size_t iPolygon = _polygons.size(); @@ -1147,9 +1186,11 @@ bool Hexahedron::compute( const Solid* solid, const IsInternalFlag intFlag ) ( nbFreeLinks < 4 && nbVertexNodes == 0 )) { // get a remaining link to start from - for ( size_t iL = 0; iL < freeLinks.size() && !curLink; ++iL ) - if (( curLink = freeLinks[ iL ] )) + for ( size_t iL = 0; iL < freeLinks.size() && !curLink; ++iL ) { + curLink = freeLinks[ iL ]; + if ( curLink ) freeLinks[ iL ] = 0; + } polygon._links.push_back( *curLink ); --nbFreeLinks; do @@ -1170,109 +1211,70 @@ bool Hexahedron::compute( const Solid* solid, const IsInternalFlag intFlag ) else // there are intersections with EDGEs { // get a remaining link to start from, one lying on minimal nb of FACEs - { - typedef pair< TGeomID, int > TFaceOfLink; - TFaceOfLink faceOfLink( -1, -1 ); - TFaceOfLink facesOfLink[3] = { faceOfLink, faceOfLink, faceOfLink }; - for ( size_t iL = 0; iL < freeLinks.size(); ++iL ) - if ( freeLinks[ iL ] ) - { - faces = freeLinks[ iL ]->GetNotUsedFace( usedFaceIDs ); - if ( faces.size() == 1 ) - { - faceOfLink = TFaceOfLink( faces[0], iL ); - if ( !freeLinks[ iL ]->HasEdgeNodes() ) - break; - facesOfLink[0] = faceOfLink; - } - else if ( facesOfLink[0].first < 0 ) - { - faceOfLink = TFaceOfLink(( faces.empty() ? -1 : faces[0]), iL ); - facesOfLink[ 1 + faces.empty() ] = faceOfLink; - } - } - for ( int i = 0; faceOfLink.first < 0 && i < 3; ++i ) - faceOfLink = facesOfLink[i]; + TFaceOfLink faceOfLink = findStartLink(freeLinks, usedFaceIDs); + curFace = faceOfLink.first; + curLink = freeLinks[ faceOfLink.second ]; + freeLinks[ faceOfLink.second ] = 0; - if ( faceOfLink.first < 0 ) // all faces used - { - for ( size_t iL = 0; iL < freeLinks.size() && faceOfLink.first < 1; ++iL ) - if (( curLink = freeLinks[ iL ])) - { - faceOfLink.first = - curLink->FirstNode()->IsLinked( curLink->LastNode()->_intPoint ); - faceOfLink.second = iL; - } - usedFaceIDs.clear(); - } - curFace = faceOfLink.first; - curLink = freeLinks[ faceOfLink.second ]; - freeLinks[ faceOfLink.second ] = 0; - } usedFaceIDs.insert( curFace ); polygon._links.push_back( *curLink ); --nbFreeLinks; // find all links lying on a curFace - do - { + do { // go forward from curLink curNode = curLink->LastNode(); curLink = 0; - for ( size_t iL = 0; iL < freeLinks.size() && !curLink; ++iL ) + for ( size_t iL = 0; iL < freeLinks.size() && !curLink; ++iL ) { if ( freeLinks[ iL ] && - freeLinks[ iL ]->FirstNode() == curNode && - freeLinks[ iL ]->LastNode()->IsOnFace( curFace )) - { + freeLinks[ iL ]->FirstNode() == curNode && + freeLinks[ iL ]->LastNode()->IsOnFace( curFace )) { curLink = freeLinks[ iL ]; freeLinks[ iL ] = 0; polygon._links.push_back( *curLink ); --nbFreeLinks; } + } } while ( curLink ); std::reverse( polygon._links.begin(), polygon._links.end() ); curLink = & polygon._links.back(); - do - { + do { // go backward from curLink curNode = curLink->FirstNode(); curLink = 0; - for ( size_t iL = 0; iL < freeLinks.size() && !curLink; ++iL ) + for ( size_t iL = 0; iL < freeLinks.size() && !curLink; ++iL ) { if ( freeLinks[ iL ] && - freeLinks[ iL ]->LastNode() == curNode && - freeLinks[ iL ]->FirstNode()->IsOnFace( curFace )) - { + freeLinks[ iL ]->LastNode() == curNode && + freeLinks[ iL ]->FirstNode()->IsOnFace( curFace )) { curLink = freeLinks[ iL ]; freeLinks[ iL ] = 0; polygon._links.push_back( *curLink ); --nbFreeLinks; } + } } while ( curLink ); curNode = polygon._links.back().FirstNode(); - if ( polygon._links[0].LastNode() != curNode ) - { - if ( nbVertexNodes > 0 ) - { + if ( polygon._links[0].LastNode() != curNode ) { + if ( nbVertexNodes > 0 ) { // add links with _vIntNodes if not already used - chainNodes.clear(); - for ( size_t iN = 0; iN < _vIntNodes.size(); ++iN ) + vector< _Node* > chainNodes; + //chainNodes.clear(); + for ( size_t iN = 0; iN < _vIntNodes.size(); ++iN ) { if ( !_vIntNodes[ iN ]->IsUsedInFace() && - _vIntNodes[ iN ]->IsOnFace( curFace )) - { + _vIntNodes[ iN ]->IsOnFace( curFace )) { _vIntNodes[ iN ]->_usedInFace = &polygon; chainNodes.push_back( _vIntNodes[ iN ] ); } + } if ( chainNodes.size() > 1 && - curFace != _grid->PseudoIntExtFaceID() ) /////// TODO - { + curFace != _grid->PseudoIntExtFaceID() ) { /////// TODO sortVertexNodes( chainNodes, curNode, curFace ); } - for ( size_t i = 0; i < chainNodes.size(); ++i ) - { + for ( size_t i = 0; i < chainNodes.size(); ++i ) { polygon.AddPolyLink( chainNodes[ i ], curNode ); curNode = chainNodes[ i ]; freeLinks.push_back( &polygon._links.back() ); @@ -1290,7 +1292,7 @@ bool Hexahedron::compute( const Solid* solid, const IsInternalFlag intFlag ) } // if there are intersections with EDGEs if ( polygon._links.size() < 2 || - polygon._links[0].LastNode() != polygon._links.back().FirstNode() ) + polygon._links[0].LastNode() != polygon._links.back().FirstNode() ) { _polygons.clear(); break; // closed polygon not found -> invalid polyhedron @@ -1311,7 +1313,7 @@ bool Hexahedron::compute( const Solid* solid, const IsInternalFlag intFlag ) if ( iPolygon == _polygons.size()-1 ) _polygons.pop_back(); } - else // polygon._links.size() >= 2 + else // polygon._links.size() > 2 { // add polygon to its links for ( size_t iL = 0; iL < polygon._links.size(); ++iL ) @@ -1322,7 +1324,7 @@ bool Hexahedron::compute( const Solid* solid, const IsInternalFlag intFlag ) if ( /*hasEdgeIntersections &&*/ iPolygon == _polygons.size() - 1 ) { // check that a polygon does not lie on a hexa side - coplanarPolyg = 0; + _Face* coplanarPolyg = 0; for ( size_t iL = 0; iL < polygon._links.size() && !coplanarPolyg; ++iL ) { if ( polygon._links[ iL ].NbFaces() < 2 ) @@ -1332,11 +1334,13 @@ bool Hexahedron::compute( const Solid* solid, const IsInternalFlag intFlag ) size_t iL2; coplanarPolyg = polygon._links[ iL ]._link->_faces[0]; for ( iL2 = iL + 1; iL2 < polygon._links.size(); ++iL2 ) + { if ( polygon._links[ iL2 ]._link->_faces[0] == coplanarPolyg && !coplanarPolyg->IsPolyLink( polygon._links[ iL ]) && !coplanarPolyg->IsPolyLink( polygon._links[ iL2 ]) && coplanarPolyg < & _polygons[ nbQuadPolygons ]) break; + } if ( iL2 == polygon._links.size() ) coplanarPolyg = 0; } @@ -1350,11 +1354,13 @@ bool Hexahedron::compute( const Solid* solid, const IsInternalFlag intFlag ) // fill freeLinks with links not shared by coplanarPolyg and polygon for ( size_t iL = 0; iL < polygon._links.size(); ++iL ) + { if ( polygon._links[ iL ]._link->_faces[1] && polygon._links[ iL ]._link->_faces[0] != coplanarPolyg ) { _Face* p = polygon._links[ iL ]._link->_faces[0]; for ( size_t iL2 = 0; iL2 < p->_links.size(); ++iL2 ) + { if ( p->_links[ iL2 ]._link == polygon._links[ iL ]._link ) { freeLinks.push_back( & p->_links[ iL2 ] ); @@ -1362,8 +1368,11 @@ bool Hexahedron::compute( const Solid* solid, const IsInternalFlag intFlag ) freeLinks.back()->RemoveFace( &polygon ); break; } + } } + } for ( size_t iL = 0; iL < coplanarPolyg->_links.size(); ++iL ) + { if ( coplanarPolyg->_links[ iL ]._link->_faces[1] && coplanarPolyg->_links[ iL ]._link->_faces[1] != &polygon ) { @@ -1371,6 +1380,7 @@ bool Hexahedron::compute( const Solid* solid, const IsInternalFlag intFlag ) if ( p == coplanarPolyg ) p = coplanarPolyg->_links[ iL ]._link->_faces[1]; for ( size_t iL2 = 0; iL2 < p->_links.size(); ++iL2 ) + { if ( p->_links[ iL2 ]._link == coplanarPolyg->_links[ iL ]._link ) { // set links of coplanarPolyg in place of used freeLinks @@ -1395,9 +1405,12 @@ bool Hexahedron::compute( const Solid* solid, const IsInternalFlag intFlag ) } break; } + } } + } // set coplanarPolyg to be re-created next for ( size_t iP = 0; iP < _polygons.size(); ++iP ) + { if ( coplanarPolyg == & _polygons[ iP ] ) { iPolygon = iP; @@ -1405,6 +1418,7 @@ bool Hexahedron::compute( const Solid* solid, const IsInternalFlag intFlag ) _polygons[ iPolygon ]._polyLinks.clear(); break; } + } _polygons.pop_back(); usedFaceIDs.erase( curFace ); continue; @@ -1415,6 +1429,7 @@ bool Hexahedron::compute( const Solid* solid, const IsInternalFlag intFlag ) } // end of case ( polygon._links.size() > 2 ) } // while ( nbFreeLinks > 0 ) + // std::cout << "3\n"; for ( auto & face_ip : tmpAddedFace ) { @@ -1435,13 +1450,14 @@ bool Hexahedron::compute( const Solid* solid, const IsInternalFlag intFlag ) _hasTooSmall = ! checkPolyhedronSize( intFlag & IS_CUT_BY_INTERNAL_FACE, volSize ); for ( size_t i = 0; i < 8; ++i ) + { if ( _hexNodes[ i ]._intPoint == &ipTmp ) _hexNodes[ i ]._intPoint = 0; + } if ( _hasTooSmall ) return false; // too small volume - // Try to find out names of no-name polygons (issue # 19887) if ( _grid->IsToRemoveExcessEntities() && _polygons.back()._name == SMESH_Block::ID_NONE ) { @@ -1474,8 +1490,7 @@ bool Hexahedron::compute( const Solid* solid, const IsInternalFlag intFlag ) } } - - /* This call is irrelevant here because _volumeDefs datas are were not filled! + /* This call is irrelevant here because _volumeDefs datas were not filled! or .... it is potentialy filled by other thread?? */ _volumeDefs._nodes.clear(); _volumeDefs._quantities.clear(); @@ -4223,4 +4238,4 @@ bool Hexahedron::_SplitIterator::Next() } } return More(); -} \ No newline at end of file +} diff --git a/src/StdMeshers/StdMeshers_Cartesian_3D_Hexahedron.hxx b/src/StdMeshers/StdMeshers_Cartesian_3D_Hexahedron.hxx index c87877b4d..ec1e1d0fe 100644 --- a/src/StdMeshers/StdMeshers_Cartesian_3D_Hexahedron.hxx +++ b/src/StdMeshers/StdMeshers_Cartesian_3D_Hexahedron.hxx @@ -35,9 +35,9 @@ #include "SMESH_StdMeshers.hxx" #include "StdMeshers_Cartesian_3D_Grid.hxx" -using namespace std; -using namespace SMESH; -namespace +namespace StdMeshers +{ +namespace Cartesian3D { // -------------------------------------------------------------------------- /*! @@ -49,14 +49,14 @@ namespace int _dInd[4][3]; size_t _nbCells[3]; int _i,_j,_k; - Grid* _grid; + StdMeshers::Cartesian3D::Grid* _grid; - CellsAroundLink( Grid* grid, int iDir ): + CellsAroundLink( StdMeshers::Cartesian3D::Grid* grid, int iDir ): _iDir( iDir ), _dInd{ {0,0,0}, {0,0,0}, {0,0,0}, {0,0,0} }, _nbCells{ grid->_coords[0].size() - 1, - grid->_coords[1].size() - 1, - grid->_coords[2].size() - 1 }, + grid->_coords[1].size() - 1, + grid->_coords[2].size() - 1 }, _grid( grid ) { const int iDirOther[3][2] = {{ 1,2 },{ 0,2 },{ 0,1 }}; @@ -85,419 +85,433 @@ namespace return true; } }; -} - -// -------------------------------------------------------------------------- -/*! - * \brief Class representing topology of the hexahedron and creating a mesh - * volume basing on analysis of hexahedron intersection with geometry - */ -class STDMESHERS_EXPORT Hexahedron -{ - // -------------------------------------------------------------------------------- - struct _Face; - struct _Link; - enum IsInternalFlag { IS_NOT_INTERNAL, IS_INTERNAL, IS_CUT_BY_INTERNAL_FACE }; - // -------------------------------------------------------------------------------- - struct _Node //!< node either at a hexahedron corner or at intersection - { - const SMDS_MeshNode* _node; // mesh node at hexahedron corner - const SMDS_MeshNode* _boundaryCornerNode; // missing mesh node due to hex truncation on the boundary - const B_IntersectPoint* _intPoint; - const _Face* _usedInFace; - char _isInternalFlags; - - _Node(const SMDS_MeshNode* n=0, const B_IntersectPoint* ip=0) - :_node(n), _intPoint(ip), _usedInFace(0), _isInternalFlags(0) {} - const SMDS_MeshNode* Node() const - { return ( _intPoint && _intPoint->_node ) ? _intPoint->_node : _node; } - const SMDS_MeshNode* BoundaryNode() const - { return _node ? _node : _boundaryCornerNode; } - const E_IntersectPoint* EdgeIntPnt() const - { return static_cast< const E_IntersectPoint* >( _intPoint ); } - const F_IntersectPoint* FaceIntPnt() const - { return static_cast< const F_IntersectPoint* >( _intPoint ); } - const vector< TGeomID >& faces() const { return _intPoint->_faceIDs; } - TGeomID face(size_t i) const { return _intPoint->_faceIDs[ i ]; } - void SetInternal( IsInternalFlag intFlag ) { _isInternalFlags |= intFlag; } - bool IsCutByInternal() const { return _isInternalFlags & IS_CUT_BY_INTERNAL_FACE; } - bool IsUsedInFace( const _Face* polygon = 0 ) - { - return polygon ? ( _usedInFace == polygon ) : bool( _usedInFace ); - } - TGeomID IsLinked( const B_IntersectPoint* other, - TGeomID avoidFace=-1 ) const // returns id of a common face - { - return _intPoint ? _intPoint->HasCommonFace( other, avoidFace ) : 0; - } - bool IsOnFace( TGeomID faceID ) const // returns true if faceID is found - { - return _intPoint ? _intPoint->IsOnFace( faceID ) : false; - } - size_t GetCommonFaces( const B_IntersectPoint * other, TGeomID* common ) const - { - return _intPoint && other ? _intPoint->GetCommonFaces( other, common ) : 0; - } - gp_Pnt Point() const - { - if ( const SMDS_MeshNode* n = Node() ) - return SMESH_NodeXYZ( n ); - if ( const E_IntersectPoint* eip = - dynamic_cast< const E_IntersectPoint* >( _intPoint )) - return eip->_point; - return gp_Pnt( 1e100, 0, 0 ); - } - TGeomID ShapeID() const - { - if ( const E_IntersectPoint* eip = dynamic_cast< const E_IntersectPoint* >( _intPoint )) - return eip->_shapeID; - return 0; - } - void Add( const E_IntersectPoint* ip ); - }; - // -------------------------------------------------------------------------------- - struct _Link // link connecting two _Node's - { - _Node* _nodes[2]; - _Face* _faces[2]; // polygons sharing a link - vector< const F_IntersectPoint* > _fIntPoints; // GridLine intersections with FACEs - vector< _Node* > _fIntNodes; // _Node's at _fIntPoints - vector< _Link > _splits; - _Link(): _faces{ 0, 0 } {} - }; - // -------------------------------------------------------------------------------- - struct _OrientedLink + // -------------------------------------------------------------------------- + /*! + * \brief Class representing topology of the hexahedron and creating a mesh + * volume basing on analysis of hexahedron intersection with geometry + */ + class STDMESHERS_EXPORT Hexahedron { - _Link* _link; - bool _reverse; - _OrientedLink( _Link* link=0, bool reverse=false ): _link(link), _reverse(reverse) {} - void Reverse() { _reverse = !_reverse; } - size_t NbResultLinks() const { return _link->_splits.size(); } - _OrientedLink ResultLink(int i) const + // -------------------------------------------------------------------------------- + struct _Face; + struct _Link; + enum IsInternalFlag { IS_NOT_INTERNAL, IS_INTERNAL, IS_CUT_BY_INTERNAL_FACE }; + // -------------------------------------------------------------------------------- + struct _Node //!< node either at a hexahedron corner or at intersection { - return _OrientedLink(&_link->_splits[_reverse ? NbResultLinks()-i-1 : i],_reverse); - } - _Node* FirstNode() const { return _link->_nodes[ _reverse ]; } - _Node* LastNode() const { return _link->_nodes[ !_reverse ]; } - operator bool() const { return _link; } - vector< TGeomID > GetNotUsedFace(const set& usedIDs ) const // returns supporting FACEs - { - vector< TGeomID > faces; - const B_IntersectPoint *ip0, *ip1; - if (( ip0 = _link->_nodes[0]->_intPoint ) && - ( ip1 = _link->_nodes[1]->_intPoint )) + const SMDS_MeshNode* _node; // mesh node at hexahedron corner + const SMDS_MeshNode* _boundaryCornerNode; // missing mesh node due to hex truncation on the boundary + const StdMeshers::Cartesian3D::B_IntersectPoint* _intPoint; + const _Face* _usedInFace; + char _isInternalFlags; + + _Node(const SMDS_MeshNode* n=0, const StdMeshers::Cartesian3D::B_IntersectPoint* ip=0) + :_node(n), _intPoint(ip), _usedInFace(0), _isInternalFlags(0) {} + const SMDS_MeshNode* Node() const + { return ( _intPoint && _intPoint->_node ) ? _intPoint->_node : _node; } + const SMDS_MeshNode* BoundaryNode() const + { return _node ? _node : _boundaryCornerNode; } + const StdMeshers::Cartesian3D::E_IntersectPoint* EdgeIntPnt() const + { return static_cast< const StdMeshers::Cartesian3D::E_IntersectPoint* >( _intPoint ); } + const StdMeshers::Cartesian3D::F_IntersectPoint* FaceIntPnt() const + { return static_cast< const StdMeshers::Cartesian3D::F_IntersectPoint* >( _intPoint ); } + const std::vector< StdMeshers::Cartesian3D::TGeomID >& faces() const { return _intPoint->_faceIDs; } + StdMeshers::Cartesian3D::TGeomID face(size_t i) const { return _intPoint->_faceIDs[ i ]; } + void SetInternal( IsInternalFlag intFlag ) { _isInternalFlags |= intFlag; } + bool IsCutByInternal() const { return _isInternalFlags & IS_CUT_BY_INTERNAL_FACE; } + bool IsUsedInFace( const _Face* polygon = 0 ) { - for ( size_t i = 0; i < ip0->_faceIDs.size(); ++i ) - if ( ip1->IsOnFace ( ip0->_faceIDs[i] ) && - !usedIDs.count( ip0->_faceIDs[i] ) ) - faces.push_back( ip0->_faceIDs[i] ); + return polygon ? ( _usedInFace == polygon ) : bool( _usedInFace ); } - return faces; - } - bool HasEdgeNodes() const - { - return ( dynamic_cast< const E_IntersectPoint* >( _link->_nodes[0]->_intPoint ) || - dynamic_cast< const E_IntersectPoint* >( _link->_nodes[1]->_intPoint )); - } - int NbFaces() const - { - return !_link->_faces[0] ? 0 : 1 + bool( _link->_faces[1] ); - } - void AddFace( _Face* f ) - { - if ( _link->_faces[0] ) + StdMeshers::Cartesian3D::TGeomID IsLinked( const StdMeshers::Cartesian3D::B_IntersectPoint* other, + StdMeshers::Cartesian3D::TGeomID avoidFace=-1 ) const // returns id of a common face { - _link->_faces[1] = f; + return _intPoint ? _intPoint->HasCommonFace( other, avoidFace ) : 0; } - else + bool IsOnFace( StdMeshers::Cartesian3D::TGeomID faceID ) const // returns true if faceID is found { - _link->_faces[0] = f; - _link->_faces[1] = 0; + return _intPoint ? _intPoint->IsOnFace( faceID ) : false; } - } - void RemoveFace( _Face* f ) + size_t GetCommonFaces( const StdMeshers::Cartesian3D::B_IntersectPoint * other, + StdMeshers::Cartesian3D::TGeomID* common ) const + { + return _intPoint && other ? _intPoint->GetCommonFaces( other, common ) : 0; + } + gp_Pnt Point() const + { + if ( const SMDS_MeshNode* n = Node() ) + return SMESH_NodeXYZ( n ); + if ( const StdMeshers::Cartesian3D::E_IntersectPoint* eip = + dynamic_cast< const StdMeshers::Cartesian3D::E_IntersectPoint* >( _intPoint )) + return eip->_point; + return gp_Pnt( 1e100, 0, 0 ); + } + StdMeshers::Cartesian3D::TGeomID ShapeID() const + { + if ( const StdMeshers::Cartesian3D::E_IntersectPoint* eip = + dynamic_cast< const StdMeshers::Cartesian3D::E_IntersectPoint* >( _intPoint )) + return eip->_shapeID; + return 0; + } + + void Add( const StdMeshers::Cartesian3D::E_IntersectPoint* ip ); + }; + + // -------------------------------------------------------------------------------- + struct _Link // link connecting two _Node's { - if ( !_link->_faces[0] ) return; + _Node* _nodes[2]; + _Face* _faces[2]; // polygons sharing a link + std::vector< const StdMeshers::Cartesian3D::F_IntersectPoint* > _fIntPoints; // GridLine intersections with FACEs + std::vector< _Node* > _fIntNodes; // _Node's at _fIntPoints + std::vector< _Link > _splits; + _Link(): _faces{ 0, 0 } {} + }; - if ( _link->_faces[1] == f ) + // -------------------------------------------------------------------------------- + struct _OrientedLink + { + _Link* _link; + bool _reverse; + _OrientedLink( _Link* link=0, bool reverse=false ): _link(link), _reverse(reverse) {} + void Reverse() { _reverse = !_reverse; } + size_t NbResultLinks() const { return _link->_splits.size(); } + _OrientedLink ResultLink(int i) const { - _link->_faces[1] = 0; + return _OrientedLink(&_link->_splits[_reverse ? NbResultLinks()-i-1 : i],_reverse); } - else if ( _link->_faces[0] == f ) + _Node* FirstNode() const { return _link->_nodes[ _reverse ]; } + _Node* LastNode() const { return _link->_nodes[ !_reverse ]; } + operator bool() const { return _link; } + + // returns supporting FACEs + std::vector< StdMeshers::Cartesian3D::TGeomID > GetNotUsedFaces + (const std::set& usedIDs ) const { - _link->_faces[0] = 0; - if ( _link->_faces[1] ) + std::vector< StdMeshers::Cartesian3D::TGeomID > faces; + const StdMeshers::Cartesian3D::B_IntersectPoint *ip0, *ip1; + if (( ip0 = _link->_nodes[0]->_intPoint ) && + ( ip1 = _link->_nodes[1]->_intPoint )) { - _link->_faces[0] = _link->_faces[1]; + for ( size_t i = 0; i < ip0->_faceIDs.size(); ++i ) + if ( ip1->IsOnFace ( ip0->_faceIDs[i] ) && + !usedIDs.count( ip0->_faceIDs[i] ) ) + faces.push_back( ip0->_faceIDs[i] ); + } + return faces; + } + + bool HasEdgeNodes() const + { + return ( dynamic_cast< const StdMeshers::Cartesian3D::E_IntersectPoint* >( _link->_nodes[0]->_intPoint ) || + dynamic_cast< const StdMeshers::Cartesian3D::E_IntersectPoint* >( _link->_nodes[1]->_intPoint )); + } + int NbFaces() const + { + return !_link->_faces[0] ? 0 : 1 + bool( _link->_faces[1] ); + } + void AddFace( _Face* f ) + { + if ( _link->_faces[0] ) { + _link->_faces[1] = f; + } + else { + _link->_faces[0] = f; _link->_faces[1] = 0; } } - } - }; - // -------------------------------------------------------------------------------- - struct _SplitIterator //! set to _hexLinks splits on one side of INTERNAL FACEs - { - struct _Split // data of a link split - { - int _linkID; // hex link ID - _Node* _nodes[2]; - int _iCheckIteration; // iteration where split is tried as Hexahedron split - _Link* _checkedSplit; // split set to hex links - bool _isUsed; // used in a volume + void RemoveFace( const _Face* f ) + { + if ( !_link->_faces[0] ) return; - _Split( _Link & split, int iLink ): - _linkID( iLink ), _nodes{ split._nodes[0], split._nodes[1] }, - _iCheckIteration( 0 ), _isUsed( false ) - {} - bool IsCheckedOrUsed( bool used ) const { return used ? _isUsed : _iCheckIteration > 0; } - }; - _Link* _hexLinks; - std::vector< _Split > _splits; - int _iterationNb; - size_t _nbChecked; - size_t _nbUsed; - std::vector< _Node* > _freeNodes; // nodes reached while composing a split set - - _SplitIterator( _Link* hexLinks ): - _hexLinks( hexLinks ), _iterationNb(0), _nbChecked(0), _nbUsed(0) - { - _freeNodes.reserve( 12 ); - _splits.reserve( 24 ); - for ( int iL = 0; iL < 12; ++iL ) - for ( size_t iS = 0; iS < _hexLinks[ iL ]._splits.size(); ++iS ) - _splits.emplace_back( _hexLinks[ iL ]._splits[ iS ], iL ); - Next(); - } - bool More() const { return _nbUsed < _splits.size(); } - bool Next(); - }; - // -------------------------------------------------------------------------------- - struct _Face - { - SMESH_Block::TShapeID _name; - vector< _OrientedLink > _links; // links on GridLine's - vector< _Link > _polyLinks; // links added to close a polygonal face - vector< _Node* > _eIntNodes; // nodes at intersection with EDGEs - - _Face():_name( SMESH_Block::ID_NONE ) - {} - bool IsPolyLink( const _OrientedLink& ol ) - { - return _polyLinks.empty() ? false : - ( &_polyLinks[0] <= ol._link && ol._link <= &_polyLinks.back() ); - } - void AddPolyLink(_Node* n0, _Node* n1, _Face* faceToFindEqual=0) - { - if ( faceToFindEqual && faceToFindEqual != this ) { - for ( size_t iL = 0; iL < faceToFindEqual->_polyLinks.size(); ++iL ) - if ( faceToFindEqual->_polyLinks[iL]._nodes[0] == n1 && - faceToFindEqual->_polyLinks[iL]._nodes[1] == n0 ) - { - _links.push_back - ( _OrientedLink( & faceToFindEqual->_polyLinks[iL], /*reverse=*/true )); - return; + if ( _link->_faces[1] == f ) { + _link->_faces[1] = 0; + } + else if ( _link->_faces[0] == f ) { + _link->_faces[0] = 0; + if ( _link->_faces[1] ) { + _link->_faces[0] = _link->_faces[1]; + _link->_faces[1] = 0; } + } } - _Link l; - l._nodes[0] = n0; - l._nodes[1] = n1; - _polyLinks.push_back( l ); - _links.push_back( _OrientedLink( &_polyLinks.back() )); - } - }; - // -------------------------------------------------------------------------------- - struct _volumeDef // holder of nodes of a volume mesh element - { - typedef void* _ptr; + }; - struct _nodeDef + // -------------------------------------------------------------------------------- + struct _SplitIterator //! set to _hexLinks splits on one side of INTERNAL FACEs { - const SMDS_MeshNode* _node; // mesh node at hexahedron corner - const B_IntersectPoint* _intPoint; - - _nodeDef(): _node(0), _intPoint(0) {} - _nodeDef( _Node* n ): _node( n->_node), _intPoint( n->_intPoint ) {} - const SMDS_MeshNode* Node() const - { return ( _intPoint && _intPoint->_node ) ? _intPoint->_node : _node; } - const E_IntersectPoint* EdgeIntPnt() const - { return static_cast< const E_IntersectPoint* >( _intPoint ); } - _ptr Ptr() const { return Node() ? (_ptr) Node() : (_ptr) EdgeIntPnt(); } - bool operator==(const _nodeDef& other ) const { return Ptr() == other.Ptr(); } + struct _Split // data of a link split + { + int _linkID; // hex link ID + _Node* _nodes[2]; + int _iCheckIteration; // iteration where split is tried as Hexahedron split + _Link* _checkedSplit; // split set to hex links + bool _isUsed; // used in a volume + + _Split( _Link & split, int iLink ): + _linkID( iLink ), _nodes{ split._nodes[0], split._nodes[1] }, + _iCheckIteration( 0 ), _isUsed( false ) + {} + bool IsCheckedOrUsed( bool used ) const { return used ? _isUsed : _iCheckIteration > 0; } + }; + _Link* _hexLinks; + std::vector< _Split > _splits; + int _iterationNb; + size_t _nbChecked; + size_t _nbUsed; + std::vector< _Node* > _freeNodes; // nodes reached while composing a split set + + _SplitIterator( _Link* hexLinks ): + _hexLinks( hexLinks ), _iterationNb(0), _nbChecked(0), _nbUsed(0) + { + _freeNodes.reserve( 12 ); + _splits.reserve( 24 ); + for ( int iL = 0; iL < 12; ++iL ) + for ( size_t iS = 0; iS < _hexLinks[ iL ]._splits.size(); ++iS ) + _splits.emplace_back( _hexLinks[ iL ]._splits[ iS ], iL ); + Next(); + } + bool More() const { return _nbUsed < _splits.size(); } + bool Next(); }; - vector< _nodeDef > _nodes; - vector< int > _quantities; - _volumeDef* _next; // to store several _volumeDefs in a chain - TGeomID _solidID; - double _size; - const SMDS_MeshElement* _volume; // new volume - std::vector _brotherVolume; // produced due to poly split - - vector< SMESH_Block::TShapeID > _names; // name of side a polygon originates from - - _volumeDef(): _next(0), _solidID(0), _size(0), _volume(0) {} - ~_volumeDef() { delete _next; } - _volumeDef( _volumeDef& other ): - _next(0), _solidID( other._solidID ), _size( other._size ), _volume( other._volume ) - { _nodes.swap( other._nodes ); _quantities.swap( other._quantities ); other._volume = 0; - _names.swap( other._names ); } - - size_t size() const { return 1 + ( _next ? _next->size() : 0 ); } // nb _volumeDef in a chain - _volumeDef* at(int index) - { return index == 0 ? this : ( _next ? _next->at(index-1) : _next ); } - - void Set( _Node** nodes, int nb ) - { _nodes.assign( nodes, nodes + nb ); } - - void SetNext( _volumeDef* vd ) - { if ( _next ) { _next->SetNext( vd ); } else { _next = vd; }} - - bool IsEmpty() const { return (( _nodes.empty() ) && - ( !_next || _next->IsEmpty() )); } - bool IsPolyhedron() const { return ( !_quantities.empty() || - ( _next && !_next->_quantities.empty() )); } - - - struct _linkDef: public std::pair<_ptr,_ptr> // to join polygons in removeExcessSideDivision() + // -------------------------------------------------------------------------------- + struct _Face { - _nodeDef _node1;//, _node2; - mutable /*const */_linkDef *_prev, *_next; - size_t _loopIndex; - - _linkDef():_prev(0), _next(0) {} + SMESH_Block::TShapeID _name; + std::vector< _OrientedLink > _links; // links on GridLine's + std::vector< _Link > _polyLinks; // links added to close a polygonal face + std::vector< _Node* > _eIntNodes; // nodes at intersection with EDGEs - void init( const _nodeDef& n1, const _nodeDef& n2, size_t iLoop ) + _Face():_name( SMESH_Block::ID_NONE ) + {} + bool IsPolyLink( const _OrientedLink& ol ) { - _node1 = n1; //_node2 = n2; - _loopIndex = iLoop; - first = n1.Ptr(); - second = n2.Ptr(); - if ( first > second ) std::swap( first, second ); + return _polyLinks.empty() ? false : + ( &_polyLinks[0] <= ol._link && ol._link <= &_polyLinks.back() ); } - void setNext( _linkDef* next ) + void AddPolyLink(_Node* n0, _Node* n1, _Face* faceToFindEqual=0) { - _next = next; - next->_prev = this; + if ( faceToFindEqual && faceToFindEqual != this ) { + for ( size_t iL = 0; iL < faceToFindEqual->_polyLinks.size(); ++iL ) + if ( faceToFindEqual->_polyLinks[iL]._nodes[0] == n1 && + faceToFindEqual->_polyLinks[iL]._nodes[1] == n0 ) + { + _links.push_back + ( _OrientedLink( & faceToFindEqual->_polyLinks[iL], /*reverse=*/true )); + return; + } + } + _Link l; + l._nodes[0] = n0; + l._nodes[1] = n1; + _polyLinks.push_back( l ); + _links.push_back( _OrientedLink( &_polyLinks.back() )); } }; - }; - // topology of a hexahedron - _Node _hexNodes [8]; - _Link _hexLinks [12]; - _Face _hexQuads [6]; - - // faces resulted from hexahedron intersection - vector< _Face > _polygons; - - // intresections with EDGEs - vector< const E_IntersectPoint* > _eIntPoints; - - // additional nodes created at intersection points - vector< _Node > _intNodes; - - // nodes inside the hexahedron (at VERTEXes) refer to _intNodes - vector< _Node* > _vIntNodes; - - // computed volume elements - _volumeDef _volumeDefs; - - Grid* _grid; - double _sideLength[3]; - int _nbCornerNodes, _nbFaceIntNodes, _nbBndNodes; - int _origNodeInd; // index of _hexNodes[0] node within the _grid - size_t _i,_j,_k; - bool _hasTooSmall; - int _cellID; - -public: - Hexahedron(Grid* grid); - int MakeElements(SMESH_MesherHelper& helper, - const map< TGeomID, vector< TGeomID > >& edge2faceIDsMap, const int numOfThreads = 1 ); - void computeElements( const Solid* solid = 0, int solidIndex = -1 ); - -private: - Hexahedron(const Hexahedron& other, size_t i, size_t j, size_t k, int cellID ); - void init( size_t i, size_t j, size_t k, const Solid* solid=0 ); - void init( size_t i ); - void setIJK( size_t i ); - /*Auxiliary methods to extract operations from monolitic compute method*/ - void defineHexahedralFaces( std::vector< _OrientedLink >& splits, std::vector<_Node*>& chainNodes, std::set< TGeomID >& concaveFaces, bool toCheckSideDivision ); - bool compute( const Solid* solid, const IsInternalFlag intFlag ); - size_t getSolids( TGeomID ids[] ); - bool isCutByInternalFace( IsInternalFlag & maxFlag ); - void addEdges(SMESH_MesherHelper& helper, - vector< Hexahedron* >& intersectedHex, - const map< TGeomID, vector< TGeomID > >& edge2faceIDsMap); - gp_Pnt findIntPoint( double u1, double proj1, double u2, double proj2, - double proj, BRepAdaptor_Curve& curve, - const gp_XYZ& axis, const gp_XYZ& origin ); - int getEntity( const E_IntersectPoint* ip, int* facets, int& sub ); - bool addIntersection( const E_IntersectPoint* ip, - vector< Hexahedron* >& hexes, - int ijk[], int dIJK[] ); - bool isQuadOnFace( const size_t iQuad ); - bool findChain( _Node* n1, _Node* n2, _Face& quad, vector<_Node*>& chainNodes ); - bool closePolygon( _Face* polygon, vector<_Node*>& chainNodes ) const; - bool findChainOnEdge( const vector< _OrientedLink >& splits, - const _OrientedLink& prevSplit, - const _OrientedLink& avoidSplit, - const std::set< TGeomID > & concaveFaces, - size_t & iS, - _Face& quad, - vector<_Node*>& chn); - int addVolumes(SMESH_MesherHelper& helper ); - void addFaces( SMESH_MesherHelper& helper, - const vector< const SMDS_MeshElement* > & boundaryVolumes ); - void addSegments( SMESH_MesherHelper& helper, - const map< TGeomID, vector< TGeomID > >& edge2faceIDsMap ); - void getVolumes( vector< const SMDS_MeshElement* > & volumes ); - void getBoundaryElems( vector< const SMDS_MeshElement* > & boundaryVolumes ); - void removeExcessSideDivision(const vector< Hexahedron* >& allHexa); - void removeExcessNodes(vector< Hexahedron* >& allHexa); - void preventVolumesOverlapping(); - TGeomID getAnyFace() const; - void cutByExtendedInternal( std::vector< Hexahedron* >& hexes, - const TColStd_MapOfInteger& intEdgeIDs ); - gp_Pnt mostDistantInternalPnt( int hexIndex, const gp_Pnt& p1, const gp_Pnt& p2 ); - bool isOutPoint( _Link& link, int iP, SMESH_MesherHelper& helper, const Solid* solid ) const; - void sortVertexNodes(vector<_Node*>& nodes, _Node* curNode, TGeomID face); - bool isInHole() const; - bool hasStrangeEdge() const; - bool checkPolyhedronSize( bool isCutByInternalFace, double & volSize ) const; - int checkPolyhedronValidity( _volumeDef* volDef, std::vector>& splitQuantities, - std::vector>& splitNodes ); - const SMDS_MeshElement* addPolyhedronToMesh( _volumeDef* volDef, SMESH_MesherHelper& helper, const std::vector& nodes, - const std::vector& quantities ); - bool addHexa (); - bool addTetra(); - bool addPenta(); - bool addPyra (); - bool debugDumpLink( _Link* link ); - _Node* findEqualNode( vector< _Node* >& nodes, - const E_IntersectPoint* ip, - const double tol2 ) - { - for ( size_t i = 0; i < nodes.size(); ++i ) - if ( nodes[i]->EdgeIntPnt() == ip || - nodes[i]->Point().SquareDistance( ip->_point ) <= tol2 ) - return nodes[i]; - return 0; - } - bool isCorner( const _Node* node ) const { return ( node >= &_hexNodes[0] && - node - &_hexNodes[0] < 8 ); } - bool hasEdgesAround( const ConcaveFace* cf ) const; - bool isImplementEdges() const { return _grid->_edgeIntPool.nbElements(); } - bool isOutParam(const double uvw[3]) const; - - typedef boost::container::flat_map< TGeomID, size_t > TID2Nb; - static void insertAndIncrement( TGeomID id, TID2Nb& id2nbMap ) - { - TID2Nb::value_type s0( id, 0 ); - TID2Nb::iterator id2nb = id2nbMap.insert( s0 ).first; - id2nb->second++; - } -}; // class Hexahedron + // -------------------------------------------------------------------------------- + struct _volumeDef // holder of nodes of a volume mesh element + { + typedef void* _ptr; + + struct _nodeDef + { + const SMDS_MeshNode* _node; // mesh node at hexahedron corner + const StdMeshers::Cartesian3D::B_IntersectPoint* _intPoint; + + _nodeDef(): _node(0), _intPoint(0) {} + _nodeDef( _Node* n ): _node( n->_node), _intPoint( n->_intPoint ) {} + const SMDS_MeshNode* Node() const + { return ( _intPoint && _intPoint->_node ) ? _intPoint->_node : _node; } + const StdMeshers::Cartesian3D::E_IntersectPoint* EdgeIntPnt() const + { return static_cast< const StdMeshers::Cartesian3D::E_IntersectPoint* >( _intPoint ); } + _ptr Ptr() const { return Node() ? (_ptr) Node() : (_ptr) EdgeIntPnt(); } + bool operator==(const _nodeDef& other ) const { return Ptr() == other.Ptr(); } + }; + + std::vector< _nodeDef > _nodes; + std::vector< int > _quantities; + _volumeDef* _next; // to store several _volumeDefs in a chain + StdMeshers::Cartesian3D::TGeomID _solidID; + double _size; + const SMDS_MeshElement* _volume; // new volume + std::vector _brotherVolume; // produced due to poly split + std::vector< SMESH_Block::TShapeID > _names; // name of side a polygon originates from + + _volumeDef(): _next(0), _solidID(0), _size(0), _volume(0) {} + ~_volumeDef() { delete _next; } + _volumeDef( _volumeDef& other ): + _next(0), _solidID( other._solidID ), _size( other._size ), _volume( other._volume ) + { _nodes.swap( other._nodes ); _quantities.swap( other._quantities ); other._volume = 0; + _names.swap( other._names ); } + + size_t size() const { return 1 + ( _next ? _next->size() : 0 ); } // nb _volumeDef in a chain + _volumeDef* at(int index) + { return index == 0 ? this : ( _next ? _next->at(index-1) : _next ); } + + void Set( _Node** nodes, int nb ) + { _nodes.assign( nodes, nodes + nb ); } + + void SetNext( _volumeDef* vd ) + { if ( _next ) { _next->SetNext( vd ); } else { _next = vd; }} + + bool IsEmpty() const { return (( _nodes.empty() ) && + ( !_next || _next->IsEmpty() )); } + bool IsPolyhedron() const { return ( !_quantities.empty() || + ( _next && !_next->_quantities.empty() )); } + + + struct _linkDef: public std::pair<_ptr,_ptr> // to join polygons in removeExcessSideDivision() + { + _nodeDef _node1;//, _node2; + mutable /*const */_linkDef *_prev, *_next; + size_t _loopIndex; + + _linkDef():_prev(0), _next(0) {} + + void init( const _nodeDef& n1, const _nodeDef& n2, size_t iLoop ) + { + _node1 = n1; //_node2 = n2; + _loopIndex = iLoop; + first = n1.Ptr(); + second = n2.Ptr(); + if ( first > second ) std::swap( first, second ); + } + void setNext( _linkDef* next ) + { + _next = next; + next->_prev = this; + } + }; + }; + + // topology of a hexahedron + _Node _hexNodes [8]; + _Link _hexLinks [12]; + _Face _hexQuads [6]; + + // faces resulted from hexahedron intersection + std::vector< _Face > _polygons; + + // intresections with EDGEs + std::vector< const StdMeshers::Cartesian3D::E_IntersectPoint* > _eIntPoints; + + // additional nodes created at intersection points + std::vector< _Node > _intNodes; + + // nodes inside the hexahedron (at VERTEXes) refer to _intNodes + std::vector< _Node* > _vIntNodes; + + // computed volume elements + _volumeDef _volumeDefs; + + StdMeshers::Cartesian3D::Grid* _grid; + double _sideLength[3]; + int _nbCornerNodes, _nbFaceIntNodes, _nbBndNodes; + int _origNodeInd; // index of _hexNodes[0] node within the _grid + size_t _i,_j,_k; + bool _hasTooSmall; + int _cellID; + + public: + Hexahedron(StdMeshers::Cartesian3D::Grid* grid); + int MakeElements(SMESH_MesherHelper& helper, + const TEdge2faceIDsMap& edge2faceIDsMap, + const int numOfThreads = 1 ); + void computeElements( const Solid* solid = 0, int solidIndex = -1 ); + + private: + Hexahedron(const Hexahedron& other, size_t i, size_t j, size_t k, int cellID ); + void init( size_t i, size_t j, size_t k, const Solid* solid=0 ); + void init( size_t i ); + void setIJK( size_t i ); + /*Auxiliary methods to extract operations from monolitic compute method*/ + void defineHexahedralFaces( const Solid* solid, const IsInternalFlag intFlag ); + bool compute( const Solid* solid, const IsInternalFlag intFlag ); + size_t getSolids( StdMeshers::Cartesian3D::TGeomID ids[] ); + bool isCutByInternalFace( IsInternalFlag & maxFlag ); + void addEdges(SMESH_MesherHelper& helper, + std::vector< Hexahedron* >& intersectedHex, + const TEdge2faceIDsMap& edge2faceIDsMap); + gp_Pnt findIntPoint( double u1, double proj1, double u2, double proj2, + double proj, BRepAdaptor_Curve& curve, + const gp_XYZ& axis, const gp_XYZ& origin ); + int getEntity( const StdMeshers::Cartesian3D::E_IntersectPoint* ip, int* facets, int& sub ); + bool addIntersection( const StdMeshers::Cartesian3D::E_IntersectPoint* ip, + std::vector< Hexahedron* >& hexes, + int ijk[], int dIJK[] ); + bool isQuadOnFace( const size_t iQuad ); + bool findChain( _Node* n1, _Node* n2, _Face& quad, std::vector<_Node*>& chainNodes ); + bool closePolygon( _Face* polygon, std::vector<_Node*>& chainNodes ) const; + bool findChainOnEdge( const std::vector< _OrientedLink >& splits, + const _OrientedLink& prevSplit, + const _OrientedLink& avoidSplit, + const std::set< StdMeshers::Cartesian3D::TGeomID > & concaveFaces, + size_t & iS, + _Face& quad, + std::vector<_Node*>& chn); + typedef std::pair< StdMeshers::Cartesian3D::TGeomID, int > TFaceOfLink; // (face, link) + static TFaceOfLink findStartLink(const std::vector< _OrientedLink* >& freeLinks, + std::set< StdMeshers::Cartesian3D::TGeomID >& usedFaceIDs); + int addVolumes( SMESH_MesherHelper& helper ); + void addFaces( SMESH_MesherHelper& helper, + const std::vector< const SMDS_MeshElement* > & boundaryVolumes ); + void addSegments( SMESH_MesherHelper& helper, + const TEdge2faceIDsMap& edge2faceIDsMap ); + void getVolumes( std::vector< const SMDS_MeshElement* > & volumes ); + void getBoundaryElems( std::vector< const SMDS_MeshElement* > & boundaryVolumes ); + void removeExcessSideDivision(const std::vector< Hexahedron* >& allHexa); + void removeExcessNodes(std::vector< Hexahedron* >& allHexa); + void preventVolumesOverlapping(); + StdMeshers::Cartesian3D::TGeomID getAnyFace() const; + void cutByExtendedInternal( std::vector< Hexahedron* >& hexes, + const TColStd_MapOfInteger& intEdgeIDs ); + gp_Pnt mostDistantInternalPnt( int hexIndex, const gp_Pnt& p1, const gp_Pnt& p2 ); + bool isOutPoint( _Link& link, int iP, SMESH_MesherHelper& helper, const Solid* solid ) const; + void sortVertexNodes(std::vector<_Node*>& nodes, + _Node* curNode, + StdMeshers::Cartesian3D::TGeomID face); + bool isInHole() const; + bool hasStrangeEdge() const; + bool checkPolyhedronSize( bool isCutByInternalFace, double & volSize ) const; + int checkPolyhedronValidity( _volumeDef* volDef, std::vector>& splitQuantities, + std::vector>& splitNodes ); + const SMDS_MeshElement* addPolyhedronToMesh( _volumeDef* volDef, + SMESH_MesherHelper& helper, + const std::vector& nodes, + const std::vector& quantities ); + bool addHexa (); + bool addTetra(); + bool addPenta(); + bool addPyra (); + bool debugDumpLink( _Link* link ); + _Node* findEqualNode( std::vector< _Node* >& nodes, + const StdMeshers::Cartesian3D::E_IntersectPoint* ip, + const double tol2 ) + { + for ( size_t i = 0; i < nodes.size(); ++i ) + if ( nodes[i]->EdgeIntPnt() == ip || + nodes[i]->Point().SquareDistance( ip->_point ) <= tol2 ) + return nodes[i]; + return 0; + } + bool isCorner( const _Node* node ) const { return ( node >= &_hexNodes[0] && + node - &_hexNodes[0] < 8 ); } + bool hasEdgesAround( const ConcaveFace* cf ) const; + bool isImplementEdges() const { return _grid->_edgeIntPool.nbElements(); } + bool isOutParam(const double uvw[3]) const; + + typedef boost::container::flat_map< StdMeshers::Cartesian3D::TGeomID, size_t > TID2Nb; + static void insertAndIncrement( StdMeshers::Cartesian3D::TGeomID id, TID2Nb& id2nbMap ) + { + TID2Nb::value_type s0( id, 0 ); + TID2Nb::iterator id2nb = id2nbMap.insert( s0 ).first; + id2nb->second++; + } + }; // class Hexahedron +} // end namespace Cartesian3D +} // end namespace StdMeshers #endif diff --git a/test/data/HexahedronTest/NRTMJ4.brep b/test/data/HexahedronTest/NRTMJ4.brep index 3ed6ccc59..7f6412ea5 100644 --- a/test/data/HexahedronTest/NRTMJ4.brep +++ b/test/data/HexahedronTest/NRTMJ4.brep @@ -1,8 +1,8 @@ DBRep_DrawableShape -CASCADE Topology V3, (c) Open Cascade +CASCADE Topology V1, (c) Matra-Datavision Locations 0 -Curve2ds 59 +Curve2ds 54 1 0 0 1 0 1 0 0 1 0 1 100 0 0 -1 @@ -11,58 +11,53 @@ Curve2ds 59 1 0 0 1 0 1 0 0 0 -1 1 0 0 0 1 -1 0 0 1 0 1 0 60 1 0 -1 0 0 0 1 1 0 0 1 0 1 16.899999999999999 50 0 1 1 0 0 1 0 -1 0 0 1 0 1 16.899999999999999 50 1 0 -1 59.899999999999999 0 0 1 1 0 0 1 0 1 76.799999999999997 50 0 1 +1 0 0 1 0 +1 0 60 1 0 +1 0 0 1 0 1 100 0 0 1 1 0 0 1 0 1 0 0 0 1 1 0 0 1 0 -1 100 0 0 1 1 0 50 1 0 -1 100 0 0 -1 +1 100 0 0 1 1 60 0 0 1 -1 0 -50 1 0 +1 100 0 0 -1 1 0 60 1 0 -1 0 0 0 1 -1 0 50 1 0 -1 16.899999999999999 50 0 1 1 0 -50 1 0 -1 0 0 1 0 -1 16.899999999999999 50 1 0 -1 59.899999999999999 0 0 1 +1 16.899999999999999 50 0 1 1 0 50 1 0 +1 16.899999999999999 50 1 0 +1 0 -50 1 0 1 76.799999999999997 50 0 1 +1 0 50 1 0 +1 0 60 1 0 +1 0 -50 1 0 1 0 0 0 1 1 0 50 1 0 -1 0 0 0 -1 1 60 0 0 1 -1 0 0 1 0 -2 0 0 1 0 -0 1 10 -2 25 25 -1 0 0 1 10 1 0 0 0 -1 -1 10 0 0 1 +2 25 25 -1 0 0 1 10 +1 0 0 1 0 1 16.899999999999999 0 0 -1 -1 0 0 0 -1 +1 10 0 0 1 1 0 0 0 1 +1 0 0 0 -1 1 59.899999999999999 0 0 -1 1 0 0 0 1 -1 59.899999999999999 0 0 -1 1 10 0 0 1 1 76.799999999999997 0 0 -1 1 0 50 1 0 2 0 0 1 0 -0 1 10 1 6.2831853071795862 -0 0 1 1 0 -0 0 1 -Curves 25 +Curves 27 1 0 0 0 0 0 1 1 0 0 100 -0 1 0 1 0 50 0 0 0 1 @@ -71,6 +66,7 @@ Curves 25 1 50 0 16.899999999999999 1 0 -0 1 50 0 16.899999999999999 0 0 1 1 50 0 76.799999999999997 1 0 -0 +1 60 0 0 0 0 1 1 0 0 100 1 0 -0 1 0 0 0 1 0 -0 1 0 50 100 1 0 -0 @@ -79,6 +75,7 @@ Curves 25 1 50 50 16.899999999999999 1 0 -0 1 50 50 16.899999999999999 0 0 1 1 50 50 76.799999999999997 1 0 -0 +1 60 50 0 0 0 1 1 0 50 0 1 0 -0 1 60 0 0 -0 1 0 2 25 25 0 0 0 -1 -1 0 -0 0 1 0 10 @@ -89,156 +86,21 @@ Curves 25 2 25 25 -50 0 0 -1 -1 0 -0 0 1 0 10 1 15.000000000000002 24.999999999999996 0 0 0 -1 Polygon3D 0 -PolygonOnTriangulations 54 -2 1 2 -p 0.6000000008 1 0 100 -2 3 4 -p 0.6000000008 1 0 100 -2 2 4 -p 0.6000000008 1 0 50 -2 1 2 -p 0.6000000008 1 0 50 -2 3 4 -p 0.6000000008 1 0 100 -2 3 4 -p 0.6000000008 1 0 100 -2 1 3 -p 0.6000000008 1 0 50 -2 1 2 -p 0.6000000008 1 0 50 -2 1 2 -p 0.6000000008 1 0 16.9 -2 1 2 -p 0.6000000008 1 0 16.9 -2 8 2 -p 0.6000000008 1 0 10 -2 1 3 -p 0.6000000008 1 0 10 -2 8 7 -p 0.6000000008 1 0 59.9 -2 1 2 -p 0.6000000008 1 0 59.9 -2 7 6 -p 0.6000000008 1 0 10 -2 1 3 -p 0.6000000008 1 0 10 -2 6 5 -p 0.6000000008 1 76.8 100 -2 1 2 -p 0.6000000008 1 76.8 100 -2 4 5 -p 0.6000000008 1 0 60 -2 1 3 -p 0.6000000008 1 0 60 -2 3 1 -p 0.6000000008 1 0 60 -2 1 3 -p 0.6000000008 1 0 60 -2 2 4 -p 0.6000000008 1 0 60 -2 4 5 -p 0.6000000008 1 0 60 -2 3 4 -p 0.6000000008 1 0 50 -2 2 4 -p 0.6000000008 1 0 50 -2 1 2 -p 0.6000000008 1 0 16.9 -2 3 4 -p 0.6000000008 1 0 16.9 -2 8 2 -p 0.6000000008 1 0 10 -2 2 4 -p 0.6000000008 1 0 10 -2 8 7 -p 0.6000000008 1 0 59.9 -2 3 4 -p 0.6000000008 1 0 59.9 -2 7 6 -p 0.6000000008 1 0 10 -2 2 4 -p 0.6000000008 1 0 10 -2 6 5 -p 0.6000000008 1 76.8 100 -2 3 4 -p 0.6000000008 1 76.8 100 -2 3 1 -p 0.6000000008 1 0 60 -2 2 4 -p 0.6000000008 1 0 60 -2 3 4 -p 0.6000000008 1 0 50 -2 1 3 -p 0.6000000008 1 0 50 -37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 -p 0.6000000008 1 0 0.174532925199433 0.349065850398866 0.523598775598299 0.698131700797732 0.872664625997165 1.0471975511966 1.22173047639603 1.39626340159546 1.5707963267949 1.74532925199433 1.91986217719376 2.0943951023932 2.26892802759263 2.44346095279206 2.61799387799149 2.79252680319093 2.96705972839036 3.14159265358979 3.31612557878922 3.49065850398866 3.66519142918809 3.83972435438752 4.01425727958696 4.18879020478639 4.36332312998582 4.53785605518525 4.71238898038469 4.88692190558412 5.06145483078355 5.23598775598299 5.41052068118242 5.58505360638185 5.75958653158128 5.93411945678072 6.10865238198015 6.28318530717959 -37 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 5 -p 0.6000000008 1 0 0.174532925199433 0.349065850398866 0.523598775598299 0.698131700797732 0.872664625997165 1.0471975511966 1.22173047639603 1.39626340159546 1.5707963267949 1.74532925199433 1.91986217719376 2.0943951023932 2.26892802759263 2.44346095279206 2.61799387799149 2.79252680319093 2.96705972839036 3.14159265358979 3.31612557878922 3.49065850398866 3.66519142918809 3.83972435438752 4.01425727958696 4.18879020478639 4.36332312998582 4.53785605518525 4.71238898038469 4.88692190558412 5.06145483078355 5.23598775598299 5.41052068118242 5.58505360638185 5.75958653158128 5.93411945678072 6.10865238198015 6.28318530717959 -2 2 4 -p 0.6000000008 1 0 50 -2 3 4 -p 0.6000000008 1 0 50 -2 1 2 -p 0.6000000008 1 0 50 -2 1 3 -p 0.6000000008 1 0 50 -2 2 4 -p 0.6000000008 1 0 50 -2 1 2 -p 0.6000000008 1 0 50 -2 3 4 -p 0.6000000008 1 0 50 -2 1 3 -p 0.6000000008 1 0 50 -37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 -p 0.6000000008 1 0 0.174532925199433 0.349065850398866 0.523598775598299 0.698131700797732 0.872664625997165 1.0471975511966 1.22173047639603 1.39626340159546 1.5707963267949 1.74532925199433 1.91986217719376 2.0943951023932 2.26892802759263 2.44346095279206 2.61799387799149 2.79252680319093 2.96705972839036 3.14159265358979 3.31612557878922 3.49065850398866 3.66519142918809 3.83972435438752 4.01425727958696 4.18879020478639 4.36332312998582 4.53785605518525 4.71238898038469 4.88692190558412 5.06145483078355 5.23598775598299 5.41052068118242 5.58505360638185 5.75958653158128 5.93411945678072 6.10865238198015 6.28318530717959 -37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 1 -p 0.6000000008 1 0 0.174532925199433 0.349065850398866 0.523598775598299 0.698131700797732 0.872664625997165 1.0471975511966 1.22173047639603 1.39626340159546 1.5707963267949 1.74532925199433 1.91986217719376 2.0943951023932 2.26892802759263 2.44346095279206 2.61799387799149 2.79252680319093 2.96705972839036 3.14159265358979 3.31612557878922 3.49065850398866 3.66519142918809 3.83972435438752 4.01425727958696 4.18879020478639 4.36332312998582 4.53785605518525 4.71238898038469 4.88692190558412 5.06145483078355 5.23598775598299 5.41052068118242 5.58505360638185 5.75958653158128 5.93411945678072 6.10865238198015 6.28318530717959 -2 38 1 -p 0.6000000008 1 0 50 -2 74 37 -p 0.6000000008 1 0 50 -Surfaces 15 +PolygonOnTriangulations 0 +Surfaces 12 1 0 0 0 1 0 -0 0 0 1 0 -1 0 1 0 0 0 -0 1 0 0 0 1 1 0 -0 1 0 0 100 0 0 1 1 0 -0 -0 1 0 1 0 50 0 -0 1 0 0 0 1 1 0 -0 1 0 0 0 0 0 1 1 0 -0 -0 1 0 1 60 0 0 1 0 -0 0 0 1 0 -1 0 -1 50 0 16.899999999999999 -0 1 0 0 0 1 1 0 -0 1 50 0 16.899999999999999 0 0 1 1 0 -0 -0 1 0 1 50 0 16.899999999999999 1 0 -0 0 0 1 0 -1 0 1 50 0 76.799999999999997 0 0 1 1 0 -0 -0 1 0 -1 50 50 16.899999999999999 -0 1 0 0 0 1 1 0 -0 +1 60 0 0 1 0 -0 0 0 1 0 -1 0 2 25 25 0 0 0 -1 -1 0 -0 0 1 0 10 -1 25 25 0 0 0 -1 -1 0 -0 0 1 0 -1 60 0 16.899999999999999 1 0 -0 0 0 1 0 -1 0 1 25 25 -50 0 0 -1 -1 0 -0 0 1 0 -Triangulations 12 -4 2 1 0 0 -0 0 0 0 0 100 0 50 0 0 50 100 0 0 100 0 0 -50 100 -50 2 1 3 2 3 4 -8 6 1 0 0 -60 0 0 60 0 16.9 0 0 0 0 0 100 60 0 100 60 0 76.8 50 0 76.8 50 0 16.9 0 60 16.9 60 0 0 100 0 100 60 76.8 60 76.8 50 16.9 50 8 1 3 2 1 8 7 3 4 7 8 3 5 6 7 5 7 4 -4 2 1 0 0 -0 0 100 0 50 100 60 0 100 60 50 100 0 0 0 50 60 0 60 50 3 2 1 4 2 3 -8 6 1 0 0 -60 50 0 60 50 16.9 0 50 0 0 50 100 60 50 100 60 50 76.8 50 50 76.8 50 50 16.9 0 60 16.9 60 0 0 100 0 100 60 76.8 60 76.8 50 16.9 50 8 1 3 2 1 8 7 3 4 7 8 3 5 6 7 5 7 4 -40 40 1 0 1.77635683940025e-15 -0 0 0 0 50 0 60 0 0 60 50 0 15 25 0 15.1519224698779 26.7364817766693 0 15.6030737921409 28.4202014332567 0 16.3397459621556 30 0 17.3395555688102 31.4278760968654 0 18.5721239031346 32.6604444311898 0 20 33.6602540378444 0 21.5797985667433 34.3969262078591 0 23.2635182233307 34.8480775301221 0 25 35 0 26.7364817766693 34.8480775301221 0 28.4202014332567 34.3969262078591 0 30 33.6602540378444 0 31.4278760968654 32.6604444311898 0 32.6604444311898 31.4278760968654 0 33.6602540378444 30 0 34.3969262078591 28.4202014332567 0 34.8480775301221 26.7364817766693 0 35 25 0 34.8480775301221 23.2635182233307 0 34.3969262078591 21.5797985667433 0 33.6602540378444 20 0 32.6604444311898 18.5721239031346 0 31.4278760968654 17.3395555688102 0 30 16.3397459621556 0 28.4202014332567 15.6030737921409 0 26.7364817766693 15.1519224698779 0 25 15 0 23.2635182233307 15.1519224698779 0 21.5797985667433 15.6030737921409 0 20 16.3397459621556 0 18.5721239031346 17.3395555688102 0 17.3395555688102 18.5721239031346 0 16.3397459621556 20 0 15.6030737921409 21.5797985667433 0 15.1519224698779 23.2635182233307 0 0 0 0 50 60 0 60 50 15 25 15.1519224698779 26.7364817766693 15.6030737921409 28.4202014332567 16.3397459621556 30 17.3395555688102 31.4278760968654 18.5721239031346 32.6604444311898 20 33.6602540378444 21.5797985667433 34.3969262078591 23.2635182233307 34.8480775301221 25 35 26.7364817766693 34.8480775301221 28.4202014332567 34.3969262078591 30 33.6602540378444 31.4278760968654 32.6604444311898 32.6604444311898 31.4278760968654 33.6602540378444 30 34.3969262078591 28.4202014332567 34.8480775301221 26.7364817766693 35 25 34.8480775301221 23.2635182233307 34.3969262078591 21.5797985667433 33.6602540378444 20 32.6604444311898 18.5721239031346 31.4278760968654 17.3395555688102 30 16.3397459621556 28.4202014332567 15.6030737921409 26.7364817766693 15.1519224698779 25 15 23.2635182233307 15.1519224698779 21.5797985667433 15.6030737921409 20 16.3397459621556 18.5721239031346 17.3395555688102 17.3395555688102 18.5721239031346 16.3397459621556 20 15.6030737921409 21.5797985667433 15.1519224698779 23.2635182233307 36 37 1 38 1 37 35 36 1 39 1 38 34 35 1 40 1 39 33 34 1 5 1 40 32 33 1 31 32 1 2 5 6 2 6 7 2 7 8 2 8 9 2 1 5 10 2 9 11 2 10 12 2 11 13 2 12 3 24 25 3 25 26 3 26 27 3 27 28 3 28 29 3 29 30 3 30 31 3 31 1 14 2 13 23 24 3 15 2 14 4 2 15 4 15 16 4 16 17 4 17 18 4 18 19 4 19 20 4 20 21 4 21 22 4 22 23 4 23 3 -4 2 1 0 0 -60 0 0 60 0 16.9 60 50 0 60 50 16.9 0 0 16.9 0 0 -50 16.9 -50 1 3 4 2 1 4 -4 2 1 0 0 -50 0 16.9 50 50 16.9 60 0 16.9 60 50 16.9 0 0 0 50 10 0 10 50 4 2 1 4 1 3 -4 2 1 0 7.105427357601e-15 -50 0 16.9 50 0 76.8 50 50 16.9 50 50 76.8 0 0 59.9 0 0 -50 59.9 -50 2 1 3 2 3 4 -4 2 1 0 0 -50 0 76.8 50 50 76.8 60 0 76.8 60 50 76.8 0 0 0 50 10 0 10 50 4 2 1 4 1 3 -4 2 1 0 0 -60 0 76.8 60 0 100 60 50 76.8 60 50 100 76.8 0 100 0 76.8 -50 100 -50 1 3 4 2 1 4 -74 72 1 0 0.0380530190825497 -15 25 -50 15.1519224698779 26.7364817766693 -50 15.6030737921409 28.4202014332567 -50 16.3397459621556 30 -50 17.3395555688102 31.4278760968654 -50 18.5721239031346 32.6604444311898 -50 20 33.6602540378444 -50 21.5797985667433 34.3969262078591 -50 23.2635182233307 34.8480775301221 -50 25 35 -50 26.7364817766693 34.8480775301221 -50 28.4202014332567 34.3969262078591 -50 30 33.6602540378444 -50 31.4278760968654 32.6604444311898 -50 32.6604444311898 31.4278760968654 -50 33.6602540378444 30 -50 34.3969262078591 28.4202014332567 -50 34.8480775301221 26.7364817766693 -50 35 25 -50 34.8480775301221 23.2635182233307 -50 34.3969262078591 21.5797985667433 -50 33.6602540378444 20 -50 32.6604444311898 18.5721239031346 -50 31.4278760968654 17.3395555688102 -50 30 16.3397459621556 -50 28.4202014332567 15.6030737921409 -50 26.7364817766693 15.1519224698779 -50 25 15 -50 23.2635182233307 15.1519224698779 -50 21.5797985667433 15.6030737921409 -50 20 16.3397459621556 -50 18.5721239031346 17.3395555688102 -50 17.3395555688102 18.5721239031346 -50 16.3397459621556 20 -50 15.6030737921409 21.5797985667433 -50 15.1519224698779 23.2635182233307 -50 15 25 -50 15 25 0 15.1519224698779 26.7364817766693 0 15.6030737921409 28.4202014332567 0 16.3397459621556 30 0 17.3395555688102 31.4278760968654 0 18.5721239031346 32.6604444311898 0 20 33.6602540378444 0 21.5797985667433 34.3969262078591 0 23.2635182233307 34.8480775301221 0 25 35 0 26.7364817766693 34.8480775301221 0 28.4202014332567 34.3969262078591 0 30 33.6602540378444 0 31.4278760968654 32.6604444311898 0 32.6604444311898 31.4278760968654 0 33.6602540378444 30 0 34.3969262078591 28.4202014332567 0 34.8480775301221 26.7364817766693 0 35 25 0 34.8480775301221 23.2635182233307 0 34.3969262078591 21.5797985667433 0 33.6602540378444 20 0 32.6604444311898 18.5721239031346 0 31.4278760968654 17.3395555688102 0 30 16.3397459621556 0 28.4202014332567 15.6030737921409 0 26.7364817766693 15.1519224698779 0 25 15 0 23.2635182233307 15.1519224698779 0 21.5797985667433 15.6030737921409 0 20 16.3397459621556 0 18.5721239031346 17.3395555688102 0 17.3395555688102 18.5721239031346 0 16.3397459621556 20 0 15.6030737921409 21.5797985667433 0 15.1519224698779 23.2635182233307 0 15 25 0 0 50 0.174532925199433 50 0.349065850398866 50 0.523598775598299 50 0.698131700797732 50 0.872664625997165 50 1.0471975511966 50 1.22173047639603 50 1.39626340159546 50 1.5707963267949 50 1.74532925199433 50 1.91986217719376 50 2.0943951023932 50 2.26892802759263 50 2.44346095279206 50 2.61799387799149 50 2.79252680319093 50 2.96705972839036 50 3.14159265358979 50 3.31612557878922 50 3.49065850398866 50 3.66519142918809 50 3.83972435438752 50 4.01425727958696 50 4.18879020478639 50 4.36332312998582 50 4.53785605518525 50 4.71238898038469 50 4.88692190558412 50 5.06145483078355 50 5.23598775598299 50 5.41052068118242 50 5.58505360638185 50 5.75958653158128 50 5.93411945678072 50 6.10865238198015 50 6.28318530717959 50 0 0 0.174532925199433 0 0.349065850398866 0 0.523598775598299 0 0.698131700797732 0 0.872664625997165 0 1.0471975511966 0 1.22173047639603 0 1.39626340159546 0 1.5707963267949 0 1.74532925199433 0 1.91986217719376 0 2.0943951023932 0 2.26892802759263 0 2.44346095279206 0 2.61799387799149 0 2.79252680319093 0 2.96705972839036 0 3.14159265358979 0 3.31612557878922 0 3.49065850398866 0 3.66519142918809 0 3.83972435438752 0 4.01425727958696 0 4.18879020478639 0 4.36332312998582 0 4.53785605518525 0 4.71238898038469 0 4.88692190558412 0 5.06145483078355 0 5.23598775598299 0 5.41052068118242 0 5.58505360638185 0 5.75958653158128 0 5.93411945678072 0 6.10865238198015 0 6.28318530717959 0 2 1 38 2 38 39 3 39 40 3 2 39 4 40 41 4 3 40 5 41 42 5 42 43 5 4 41 6 43 44 6 5 43 7 6 44 8 7 44 8 44 45 9 8 45 9 45 46 10 9 46 10 46 47 11 10 47 11 47 48 12 11 48 12 48 49 13 12 49 13 49 50 14 13 50 14 50 51 15 14 51 15 51 52 16 15 52 16 52 53 17 16 53 17 53 54 18 17 54 18 54 55 18 55 56 19 18 56 20 56 57 20 19 56 21 57 58 21 20 57 22 58 59 22 21 58 23 59 60 23 22 59 24 60 61 24 23 60 25 61 62 25 24 61 26 62 63 26 25 62 27 63 64 27 64 65 27 26 63 28 27 65 29 65 66 29 28 65 30 66 67 30 29 66 31 67 68 31 30 67 32 31 68 32 68 69 33 32 69 33 69 70 34 33 70 34 70 71 35 34 71 35 71 72 36 35 72 36 72 73 37 36 73 37 73 74 -36 34 1 0 7.94410929039127e-15 -15 25 -50 15.1519224698779 26.7364817766693 -50 15.6030737921409 28.4202014332567 -50 16.3397459621556 30 -50 17.3395555688102 31.4278760968654 -50 18.5721239031346 32.6604444311898 -50 20 33.6602540378444 -50 21.5797985667433 34.3969262078591 -50 23.2635182233307 34.8480775301221 -50 25 35 -50 26.7364817766693 34.8480775301221 -50 28.4202014332567 34.3969262078591 -50 30 33.6602540378444 -50 31.4278760968654 32.6604444311898 -50 32.6604444311898 31.4278760968654 -50 33.6602540378444 30 -50 34.3969262078591 28.4202014332567 -50 34.8480775301221 26.7364817766693 -50 35 25 -50 34.8480775301221 23.2635182233307 -50 34.3969262078591 21.5797985667433 -50 33.6602540378444 20 -50 32.6604444311898 18.5721239031346 -50 31.4278760968654 17.3395555688102 -50 30 16.3397459621556 -50 28.4202014332567 15.6030737921409 -50 26.7364817766693 15.1519224698779 -50 25 15 -50 23.2635182233307 15.1519224698779 -50 21.5797985667433 15.6030737921409 -50 20 16.3397459621556 -50 18.5721239031346 17.3395555688102 -50 17.3395555688102 18.5721239031346 -50 16.3397459621556 20 -50 15.6030737921409 21.5797985667433 -50 15.1519224698779 23.2635182233307 -50 10 -1.77635683940025e-15 9.84807753012208 1.7364817766693 9.39692620785909 3.42020143325669 8.66025403784439 5 7.66044443118978 6.42787609686539 6.42787609686539 7.66044443118978 5 8.66025403784439 3.42020143325669 9.39692620785908 1.7364817766693 9.84807753012208 0 10 -1.7364817766693 9.84807753012208 -3.42020143325669 9.39692620785909 -5 8.66025403784439 -6.42787609686539 7.66044443118978 -7.66044443118978 6.4278760968654 -8.66025403784438 5.00000000000001 -9.39692620785908 3.4202014332567 -9.84807753012208 1.73648177666932 -10 1.4210854715202e-14 -9.84807753012208 -1.73648177666929 -9.39692620785909 -3.42020143325667 -8.6602540378444 -4.99999999999998 -7.6604444311898 -6.42787609686538 -6.42787609686541 -7.66044443118977 -5.00000000000002 -8.66025403784437 -3.42020143325671 -9.39692620785907 -1.73648177666933 -9.84807753012208 -2.8421709430404e-14 -10 1.73648177666927 -9.84807753012209 3.42020143325666 -9.3969262078591 4.99999999999997 -8.6602540378444 6.42787609686536 -7.6604444311898 7.66044443118976 -6.42787609686542 8.66025403784437 -5.00000000000004 9.39692620785907 -3.42020143325673 9.84807753012207 -1.73648177666935 22 23 24 20 21 22 28 26 27 18 19 20 18 24 25 18 22 24 18 20 22 30 28 29 31 25 26 31 26 28 31 28 30 33 31 32 13 14 15 13 15 16 13 16 17 13 17 18 35 33 34 35 31 33 36 18 25 36 25 31 36 31 35 11 12 13 11 13 18 9 10 11 2 36 1 8 9 11 8 11 18 3 36 2 4 36 3 4 7 8 4 18 36 4 8 18 6 7 4 6 4 5 +Triangulations 0 TShapes 72 Ve @@ -260,8 +122,6 @@ Ed 1 1 0 0 100 2 1 1 0 0 100 2 2 2 0 0 100 -6 1 1 0 -6 2 2 0 0 0101000 @@ -278,8 +138,6 @@ Ed 1 2 0 0 50 2 3 1 0 0 50 2 4 3 0 0 50 -6 3 1 0 -6 4 3 0 0 0101000 @@ -296,8 +154,6 @@ Ed 1 3 0 0 100 2 5 1 0 0 100 2 6 4 0 0 100 -6 5 1 0 -6 6 4 0 0 0101000 @@ -307,8 +163,6 @@ Ed 1 4 0 0 50 2 7 1 0 0 50 2 8 5 0 0 50 -6 7 1 0 -6 8 5 0 0 0101000 @@ -319,7 +173,7 @@ Wi -70 0 -68 0 +66 0 +65 0 * Fa 0 1e-07 1 0 -2 1 + 0101000 +64 0 * Ve @@ -339,10 +193,8 @@ Ve Ed 1e-07 1 1 0 1 5 0 0 16.9 -2 9 6 0 0 16.9 -2 10 2 0 0 16.9 -6 9 2 0 -6 10 6 0 +2 9 2 0 0 16.9 +2 10 6 0 0 16.9 0 0101000 @@ -357,11 +209,8 @@ Ve Ed 1e-07 1 1 0 1 6 0 0 10 -2 11 7 0 0 10 -2 12 8 0 0 10 -2 13 2 0 0 10 -6 11 2 0 -6 12 7 0 +2 11 2 0 0 10 +2 12 7 0 0 10 0 0101000 @@ -376,11 +225,8 @@ Ve Ed 1e-07 1 1 0 1 7 0 0 59.9 -2 14 9 0 0 59.9 -2 15 7 0 0 59.9 -2 16 2 0 0 59.9 -6 13 2 0 -6 14 8 0 +2 13 2 0 0 59.9 +2 14 8 0 0 59.9 0 0101000 @@ -395,11 +241,8 @@ Ve Ed 1e-07 1 1 0 1 8 0 0 10 -2 17 7 0 0 10 -2 18 10 0 0 10 -2 19 2 0 0 10 -6 15 2 0 -6 16 9 0 +2 15 2 0 0 10 +2 16 9 0 0 10 0 0101000 @@ -413,33 +256,27 @@ Ve * Ed 1e-07 1 1 0 -1 5 0 76.8 100 -2 9 6 0 76.8 100 -2 10 2 0 76.8 100 -6 17 2 0 -6 18 10 0 +1 9 0 76.8 100 +2 17 2 0 76.8 100 +2 18 10 0 76.8 100 0 0101000 +55 0 -53 0 * Ed 1e-07 1 1 0 -1 9 0 0 60 -2 20 2 0 0 60 -2 21 3 0 0 60 -6 19 2 0 -6 20 3 0 +1 10 0 0 60 +2 19 2 0 0 60 +2 20 3 0 0 60 0 0101000 -53 0 +72 0 * Ed 1e-07 1 1 0 -1 10 0 0 60 -2 22 2 0 0 60 -2 23 5 0 0 60 -6 21 2 0 -6 22 5 0 +1 11 0 0 60 +2 21 2 0 0 60 +2 22 5 0 0 60 0 0101000 @@ -450,7 +287,7 @@ Wi -60 0 +58 0 -56 0 -54 0 -52 0 +51 0 +70 0 -50 0 * Fa 0 1e-07 2 0 -2 2 + 0101000 +49 0 * Ve @@ -462,22 +299,18 @@ Ve * Ed 1e-07 1 1 0 -1 11 0 0 60 +1 12 0 0 60 +2 23 3 0 0 60 2 24 4 0 0 60 -2 25 3 0 0 60 -6 23 3 0 -6 24 4 0 0 0101000 -47 0 +69 0 * Ed 1e-07 1 1 0 -1 12 0 0 50 -2 26 6 0 0 50 -2 27 3 0 0 50 -6 25 3 0 -6 26 10 0 +1 13 0 0 50 +2 25 3 0 0 50 +2 26 10 0 0 50 0 0101000 @@ -488,7 +321,7 @@ Wi -68 0 -46 0 +45 0 +51 0 * Fa 0 1e-07 3 0 -2 3 + 0101000 +44 0 * Ve @@ -507,11 +340,9 @@ Ve * Ed 1e-07 1 1 0 -1 13 0 0 16.9 +1 14 0 0 16.9 +2 27 4 0 0 16.9 2 28 6 0 0 16.9 -2 29 4 0 0 16.9 -6 27 4 0 -6 28 6 0 0 0101000 @@ -525,12 +356,9 @@ Ve * Ed 1e-07 1 1 0 -1 14 0 0 10 -2 30 11 0 0 10 -2 31 8 0 0 10 -2 32 4 0 0 10 -6 29 4 0 -6 30 7 0 +1 15 0 0 10 +2 29 4 0 0 10 +2 30 7 0 0 10 0 0101000 @@ -544,12 +372,9 @@ Ve * Ed 1e-07 1 1 0 -1 15 0 0 59.9 -2 33 9 0 0 59.9 -2 34 11 0 0 59.9 -2 35 4 0 0 59.9 -6 31 4 0 -6 32 8 0 +1 16 0 0 59.9 +2 31 4 0 0 59.9 +2 32 8 0 0 59.9 0 0101000 @@ -563,34 +388,27 @@ Ve * Ed 1e-07 1 1 0 -1 16 0 0 10 -2 36 11 0 0 10 -2 37 10 0 0 10 -2 38 4 0 0 10 -6 33 4 0 -6 34 9 0 +1 17 0 0 10 +2 33 4 0 0 10 +2 34 9 0 0 10 0 0101000 -35 0 +37 0 * Ed 1e-07 1 1 0 -1 13 0 76.8 100 -2 28 6 0 76.8 100 -2 29 4 0 76.8 100 -6 35 4 0 -6 36 10 0 +1 18 0 76.8 100 +2 35 4 0 76.8 100 +2 36 10 0 76.8 100 0 0101000 +35 0 -47 0 * Ed 1e-07 1 1 0 -1 17 0 0 60 -2 39 4 0 0 60 -2 40 5 0 0 60 -6 37 4 0 -6 38 5 0 +1 19 0 0 60 +2 37 4 0 0 60 +2 38 5 0 0 60 0 0101000 @@ -601,16 +419,14 @@ Wi -40 0 +38 0 -36 0 -34 0 -33 0 +46 0 +66 0 -32 0 * Fa 0 1e-07 4 0 -2 4 + 0101000 +31 0 * Ed 1e-07 1 1 0 -1 18 0 0 50 -2 41 6 0 0 50 -2 42 5 0 0 50 -6 39 5 0 -6 40 6 0 +1 20 0 0 50 +2 39 5 0 0 50 +2 40 6 0 0 50 0 0101000 @@ -628,12 +444,9 @@ Ve * Ed 1e-07 1 1 0 -1 19 0 0 6.28318530717959 -2 43 12 0 0 6.28318530717959 -2 44 13 0 0 6.28318530717959 -2 45 5 0 0 6.28318530717959 -6 41 11 0 -6 42 5 0 +1 21 0 0 6.28318530717959 +2 41 5 0 0 6.28318530717959 +2 42 11 0 0 6.28318530717959 0 0101000 @@ -644,17 +457,14 @@ Wi +26 0 * Fa 0 1e-07 5 0 -2 5 + 0101000 +28 0 +25 0 * Ed 1e-07 1 1 0 -1 20 0 0 50 -2 46 14 0 0 50 -2 47 8 0 0 50 -2 48 6 0 0 50 -6 43 6 0 -6 44 7 0 +1 22 0 0 50 +2 43 6 0 0 50 +2 44 7 0 0 50 0 0101000 @@ -665,16 +475,14 @@ Wi -60 0 -23 0 +40 0 +29 0 * Fa 0 1e-07 6 0 -2 6 + 0101000 +22 0 * Ed 1e-07 1 1 0 -1 21 0 0 50 -2 49 9 0 0 50 -2 50 8 0 0 50 -6 45 7 0 -6 46 8 0 +1 23 0 0 50 +2 45 7 0 0 50 +2 46 8 0 0 50 0 0101000 @@ -684,17 +492,15 @@ Wi 0101100 -20 0 -38 0 +23 0 +58 0 * Fa -0 1e-07 8 0 -2 7 +0 1e-07 7 0 + 0101000 +19 0 * Ed 1e-07 1 1 0 -1 22 0 0 50 -2 51 9 0 0 50 -2 52 10 0 0 50 -6 47 8 0 -6 48 9 0 +1 24 0 0 50 +2 47 8 0 0 50 +2 48 9 0 0 50 0 0101000 @@ -704,18 +510,15 @@ Wi 0101100 -56 0 -17 0 +36 0 +20 0 * Fa -0 1e-07 9 0 -2 8 +0 1e-07 8 0 + 0101000 +16 0 * Ed 1e-07 1 1 0 -1 23 0 0 50 -2 53 14 0 0 50 -2 54 10 0 0 50 -2 55 6 0 0 50 -6 49 9 0 -6 50 10 0 +1 25 0 0 50 +2 49 9 0 0 50 +2 50 10 0 0 50 0 0101000 @@ -725,8 +528,8 @@ Wi 0101100 -17 0 -34 0 +14 0 +54 0 * Fa -0 1e-07 10 0 -2 9 +0 1e-07 9 0 + 0101000 +13 0 * Wi @@ -734,8 +537,8 @@ Wi 0101100 -52 0 -45 0 +33 0 +14 0 * Fa -0 1e-07 6 0 -2 10 +0 1e-07 10 0 + 0101000 +11 0 * Ve @@ -747,20 +550,17 @@ Ve * Ed 1e-07 1 1 0 -1 24 0 0 6.28318530717959 -2 56 12 0 0 6.28318530717959 -2 57 15 0 0 6.28318530717959 -6 51 11 0 -6 52 12 0 +1 26 0 0 6.28318530717959 +2 51 11 0 0 6.28318530717959 +2 52 12 0 0 6.28318530717959 0 0101000 +9 0 -9 0 * Ed 1e-07 1 1 0 -1 25 0 0 50 -3 58 59CN 12 0 0 50 -7 53 54 11 0 +1 27 0 0 50 +3 53 54CN 11 0 0 50 0 0101000 @@ -770,8 +570,8 @@ Wi 0101100 -8 0 +7 0 +26 0 -7 0 * Fa -0 1e-07 12 0 -2 11 +0 1e-07 11 0 + 0101000 +6 0 * Wi @@ -779,8 +579,8 @@ Wi 0101100 +8 0 * Fa -0 1e-07 15 0 -2 12 +0 1e-07 12 0 + 0101000 +4 0 * Sh @@ -790,7 +590,9 @@ Sh +5 0 +3 0 * So -0100000 +1100000 +2 0 * -+1 0 \ No newline at end of file ++1 0 +0 + -- 2.39.2