0029715: Mesh - Estimate the grid size of the acceleration structure by the complexit...

The circles acceleration structure sizes are estimated using the required deflection and the complexity of the face.

diff --git a/src/BRepMesh/BRepMesh_Delaun.cxx b/src/BRepMesh/BRepMesh_Delaun.cxx
index 62ff759122532aaa76dadb6f5dc5cf42bbe872a6..a6aec38ecb0985b7645adec999d5436e1839e76c 100644 (file)
--- a/src/BRepMesh/BRepMesh_Delaun.cxx
+++ b/src/BRepMesh/BRepMesh_Delaun.cxx
@@ -106,7 +106,9 @@ BRepMesh_Delaun::BRepMesh_Delaun(
   myCircles ( theVertices.Length(), theOldMesh->Allocator() )
 {
   if ( theVertices.Length() > 2 )
+  {
     Init( theVertices );
+  }
 }
 
 //=======================================================================
@@ -119,16 +121,21 @@ BRepMesh_Delaun::BRepMesh_Delaun(
 : myMeshData( theOldMesh ),
   myCircles ( theVertexIndices.Length(), theOldMesh->Allocator() )
 {
-  if ( theVertexIndices.Length() > 2 )
-  {
-    Bnd_Box2d aBox;
-    Standard_Integer anIndex = theVertexIndices.Lower();
-    Standard_Integer anUpper = theVertexIndices.Upper();
-    for ( ; anIndex <= anUpper; ++anIndex )
-      aBox.Add( gp_Pnt2d( GetVertex( theVertexIndices( anIndex) ).Coord() ) );
+  perform( theVertexIndices );
+}
 
-    perform( aBox, theVertexIndices );
-  }
+//=======================================================================
+//function : BRepMesh_Delaun
+//purpose  : Creates the triangulation with and existent Mesh data structure
+//=======================================================================
+BRepMesh_Delaun::BRepMesh_Delaun (const Handle (BRepMesh_DataStructureOfDelaun)& theOldMesh,
+                                  BRepMesh::Array1OfInteger&                     theVertexIndices,
+                                  const Standard_Integer                         theCellsCountU,
+                                  const Standard_Integer                         theCellsCountV)
+  : myMeshData (theOldMesh),
+  myCircles (theVertexIndices.Length (), theOldMesh->Allocator ())
+{
+  perform (theVertexIndices, theCellsCountU, theCellsCountV);
 }
 
 //=======================================================================
@@ -137,7 +144,6 @@ BRepMesh_Delaun::BRepMesh_Delaun(
 //=======================================================================
 void BRepMesh_Delaun::Init(BRepMesh::Array1OfVertexOfDelaun& theVertices)
 {
-  Bnd_Box2d aBox;
   Standard_Integer aLowerIdx  = theVertices.Lower();
   Standard_Integer anUpperIdx = theVertices.Upper();
   BRepMesh::Array1OfInteger aVertexIndexes( aLowerIdx, anUpperIdx );
@@ -145,35 +151,62 @@ void BRepMesh_Delaun::Init(BRepMesh::Array1OfVertexOfDelaun& theVertices)
   Standard_Integer anIndex = aLowerIdx;
   for ( ; anIndex <= anUpperIdx; ++anIndex )
   {
-    aBox.Add( gp_Pnt2d( theVertices( anIndex ).Coord() ) );
     aVertexIndexes( anIndex ) = myMeshData->AddNode( theVertices( anIndex ) );
   }
 
-  perform( aBox, aVertexIndexes );
+  perform( aVertexIndexes );
 }
 
 //=======================================================================
 //function : perform
 //purpose  : Create super mesh and run triangulation procedure
 //=======================================================================
-void BRepMesh_Delaun::perform(Bnd_Box2d&                 theBndBox,
-                              BRepMesh::Array1OfInteger& theVertexIndexes)
+void BRepMesh_Delaun::perform (BRepMesh::Array1OfInteger& theVertexIndices,
+                               const Standard_Integer     theCellsCountU /* = -1 */,
+                               const Standard_Integer     theCellsCountV /* = -1 */)
 {
-  theBndBox.Enlarge( Precision );
-  superMesh( theBndBox );
+  if (theVertexIndices.Length () <= 2)
+  {
+    return;
+  }
+
+  Bnd_Box2d aBox;
+  Standard_Integer anIndex = theVertexIndices.Lower ();
+  Standard_Integer anUpper = theVertexIndices.Upper ();
+  for (; anIndex <= anUpper; ++anIndex)
+  {
+    aBox.Add (gp_Pnt2d (GetVertex (theVertexIndices (anIndex)).Coord ()));
+  }
+
+  aBox.Enlarge (Precision);
+
+  Standard_Integer aScaler = 2;
+  if ( myMeshData->NbNodes() > 100 )
+  {
+    aScaler = 5;
+  }
+  else if( myMeshData->NbNodes() > 1000 )
+  {
+    aScaler = 7;
+  }
+
+  superMesh (aBox, Max (theCellsCountU, aScaler),
+                   Max (theCellsCountV, aScaler));
 
   ComparatorOfIndexedVertexOfDelaun aCmp(myMeshData);
-  std::make_heap(theVertexIndexes.begin(), theVertexIndexes.end(), aCmp);
-  std::sort_heap(theVertexIndexes.begin(), theVertexIndexes.end(), aCmp);
+  std::make_heap(theVertexIndices.begin(), theVertexIndices.end(), aCmp);
+  std::sort_heap(theVertexIndices.begin(), theVertexIndices.end(), aCmp);
 
-  compute( theVertexIndexes );
+  compute( theVertexIndices );
 }
 
 //=======================================================================
 //function : superMesh
 //purpose  : Build the super mesh
 //=======================================================================
-void BRepMesh_Delaun::superMesh( const Bnd_Box2d& theBox )
+void BRepMesh_Delaun::superMesh (const Bnd_Box2d&       theBox,
+                                 const Standard_Integer theCellsCountU,
+                                 const Standard_Integer theCellsCountV)
 {
   Standard_Real aMinX, aMinY, aMaxX, aMaxY;
   theBox.Get( aMinX, aMinY, aMaxX, aMaxY );
@@ -185,15 +218,7 @@ void BRepMesh_Delaun::superMesh( const Bnd_Box2d& theBox )
   Standard_Real aDelta    = aDeltaX + aDeltaY;
   
   myCircles.SetMinMaxSize( gp_XY( aMinX, aMinY ), gp_XY( aMaxX, aMaxY ) );
-
-  Standard_Integer aScaler = 2;
-  if ( myMeshData->NbNodes() > 100 )
-    aScaler = 5;
-  else if( myMeshData->NbNodes() > 1000 )
-    aScaler = 7;
-
-  myCircles.SetCellSize( aDeltaX / aScaler,
-                         aDeltaY / aScaler );
+  myCircles.SetCellSize( aDeltaX / theCellsCountU, aDeltaY / theCellsCountV);
 
   mySupVert[0] = myMeshData->AddNode(
     BRepMesh_Vertex( ( aMinX + aMaxX ) / 2, aMaxY + aDeltaMax, BRepMesh_Free ) );

diff --git a/src/BRepMesh/BRepMesh_Delaun.hxx b/src/BRepMesh/BRepMesh_Delaun.hxx
index e8fcba595e5a975fa2c4b92627460fc861fc1691..8741d8b074bb46db2950316106019f1ab90ffcba 100755 (executable)
--- a/src/BRepMesh/BRepMesh_Delaun.hxx
+++ b/src/BRepMesh/BRepMesh_Delaun.hxx
@@ -52,6 +52,12 @@ public:
   Standard_EXPORT BRepMesh_Delaun (const Handle(BRepMesh_DataStructureOfDelaun)& theOldMesh,
                                    BRepMesh::Array1OfInteger&                    theVertexIndices);
 
+  //! Creates the triangulation with an existant Mesh data structure.
+  Standard_EXPORT BRepMesh_Delaun (const Handle (BRepMesh_DataStructureOfDelaun)& theOldMesh,
+                                   BRepMesh::Array1OfInteger&                     theVertexIndices,
+                                   const Standard_Integer                         theCellsCountU,
+                                   const Standard_Integer                         theCellsCountV);
+
   //! Initializes the triangulation with an array of vertices.
   Standard_EXPORT void Init (BRepMesh::Array1OfVertexOfDelaun& theVertices);
 
@@ -144,12 +150,15 @@ private:
   //! that have number of connected elements less or equal 1.
   BRepMesh::HMapOfInteger getEdgesByType (const BRepMesh_DegreeOfFreedom theEdgeType) const;
 
-  //! Create super mesh and run triangulation procedure.
-  void perform (Bnd_Box2d&                 theBndBox,
-                BRepMesh::Array1OfInteger& theVertexIndices);
+  //! Run triangulation procedure.
+  void perform (BRepMesh::Array1OfInteger& theVertexIndices,
+                const Standard_Integer     theCellsCountU = -1,
+                const Standard_Integer     theCellsCountV = -1);
 
   //! Build the super mesh.
-  void superMesh (const Bnd_Box2d& theBox);
+  void superMesh (const Bnd_Box2d&       theBox,
+                  const Standard_Integer theCellsCountU,
+                  const Standard_Integer theCellsCountV);
 
   //! Computes the triangulation and adds the vertices,
   //! edges and triangles to the Mesh data structure.

diff --git a/src/BRepMesh/BRepMesh_FastDiscretFace.cxx b/src/BRepMesh/BRepMesh_FastDiscretFace.cxx
index 89cbc895402456441cc2ee1f8fa9c0039e21c76b..f50bc58922798893871c8cb87c0cbee8dc632722 100644 (file)
--- a/src/BRepMesh/BRepMesh_FastDiscretFace.cxx
+++ b/src/BRepMesh/BRepMesh_FastDiscretFace.cxx
@@ -23,6 +23,7 @@
 #include <BRepAdaptor_HSurface.hxx>
 #include <BRepAdaptor_Curve.hxx>
 #include <Adaptor3d_IsoCurve.hxx>
+#include <Adaptor3d_HCurve.hxx>
 
 #include <BRep_ListIteratorOfListOfPointRepresentation.hxx>
 #include <BRep_PointRepresentation.hxx>
@@ -162,6 +163,119 @@ namespace
   private:
     const TopoDS_Vertex& myVertex;
   };
+
+  void ComputeErrFactors (const Standard_Real               theDeflection, 
+                          const Handle(Adaptor3d_HSurface)& theFace,
+                          Standard_Real&                    theErrFactorU,
+                          Standard_Real&                    theErrFactorV)
+  {
+    theErrFactorU = theDeflection * 10.;
+    theErrFactorV = theDeflection * 10.;
+
+    switch (theFace->GetType ())
+    {
+    case GeomAbs_Cylinder:
+    case GeomAbs_Cone:
+    case GeomAbs_Sphere:
+    case GeomAbs_Torus:
+      break;
+
+    case GeomAbs_SurfaceOfExtrusion:
+    case GeomAbs_SurfaceOfRevolution:
+      {
+        Handle (Adaptor3d_HCurve) aCurve = theFace->BasisCurve ();
+        if (aCurve->GetType () == GeomAbs_BSplineCurve && aCurve->Degree () > 2)
+        {
+          theErrFactorV /= (aCurve->Degree () * aCurve->NbKnots ());
+        }
+        break;
+      }
+    case GeomAbs_BezierSurface:
+      {
+        if (theFace->UDegree () > 2)
+        {
+          theErrFactorU /= (theFace->UDegree ());
+        }
+        if (theFace->VDegree () > 2)
+        {
+          theErrFactorV /= (theFace->VDegree ());
+        }
+        break;
+      }
+    case GeomAbs_BSplineSurface:
+      {
+        if (theFace->UDegree () > 2)
+        {
+          theErrFactorU /= (theFace->UDegree () * theFace->NbUKnots ());
+        }
+        if (theFace->VDegree () > 2)
+        {
+          theErrFactorV /= (theFace->VDegree () *  theFace->NbVKnots ());
+        }
+        break;
+      }
+
+    case GeomAbs_Plane:
+    default:
+      theErrFactorU = theErrFactorV = 1.;
+    }
+  }
+
+  void AdjustCellsCounts (const Handle(Adaptor3d_HSurface)& theFace,
+                          const Standard_Integer            theNbVertices,
+                          Standard_Integer&                 theCellsCountU,
+                          Standard_Integer&                 theCellsCountV)
+  {
+    const GeomAbs_SurfaceType aType = theFace->GetType ();
+    if (aType == GeomAbs_OtherSurface)
+    {
+      // fallback to the default behavior
+      theCellsCountU = theCellsCountV = -1;
+      return;
+    }
+
+    Standard_Real aSqNbVert = theNbVertices;
+    if (aType == GeomAbs_Plane)
+    {
+      theCellsCountU = theCellsCountV = (Standard_Integer)Ceiling (Pow (2, Log10 (aSqNbVert)));
+    }
+    else if (aType == GeomAbs_Cylinder || aType == GeomAbs_Cone)
+    {
+      theCellsCountV = (Standard_Integer)Ceiling (Pow (2, Log10 (aSqNbVert)));
+    }
+    else if (aType == GeomAbs_SurfaceOfExtrusion || aType == GeomAbs_SurfaceOfRevolution)
+    {
+      Handle (Adaptor3d_HCurve) aCurve = theFace->BasisCurve ();
+      if (aCurve->GetType () == GeomAbs_Line ||
+          (aCurve->GetType () == GeomAbs_BSplineCurve && aCurve->Degree () < 2))
+      {
+        // planar, cylindrical, conical cases
+        if (aType == GeomAbs_SurfaceOfExtrusion)
+          theCellsCountU = (Standard_Integer)Ceiling (Pow (2, Log10 (aSqNbVert)));
+        else
+          theCellsCountV = (Standard_Integer)Ceiling (Pow (2, Log10 (aSqNbVert)));
+      }
+      if (aType == GeomAbs_SurfaceOfExtrusion)
+      {
+        // V is always a line
+        theCellsCountV = (Standard_Integer)Ceiling (Pow (2, Log10 (aSqNbVert)));
+      }
+    }
+    else if (aType == GeomAbs_BezierSurface || aType == GeomAbs_BSplineSurface)
+    {
+      if (theFace->UDegree () < 2)
+      {
+        theCellsCountU = (Standard_Integer)Ceiling (Pow (2, Log10 (aSqNbVert)));
+      }
+      if (theFace->VDegree () < 2)
+      {
+        theCellsCountV = (Standard_Integer)Ceiling (Pow (2, Log10 (aSqNbVert)));
+      }
+    }
+
+    theCellsCountU = Max (theCellsCountU, 2);
+    theCellsCountV = Max (theCellsCountV, 2);
+  }
 }
 
 
@@ -359,7 +473,56 @@ void BRepMesh_FastDiscretFace::add(const Handle(BRepMesh_FaceAttribute)& theAttr
   for ( Standard_Integer i = 1; i <= nbVertices; ++i )
     tabvert_corr(i) = i;
 
-  BRepMesh_Delaun trigu(myStructure, tabvert_corr);
+  Handle (Adaptor3d_HSurface) aSurface (myAttribute->Surface ());
+  GeomAbs_SurfaceType aType = aSurface->GetType ();
+
+  while (aType == GeomAbs_OffsetSurface)
+  {
+    aSurface = aSurface->BasisSurface ();
+    aType = aSurface->GetType ();
+  }
+
+  // For better meshing performance we try to estimate the acceleration circles grid structure sizes:
+  // For each parametric direction (U, V) we estimate firstly an approximate distance between the future points -
+  // this estimation takes into account the required face deflection and the complexity of the face.
+  // Particularly, the complexity of the faces based on BSpline curves and surfaces requires much more points.
+  // At the same time, for planar faces and linear parts of the arbitrary surfaces usually no intermediate points
+  // are necessary.
+  // The general idea for each parametric direction:
+  // cells_count = 2 ^ log10 ( estimated_points_count )
+  // For linear parametric direction we fall back to the initial vertex count:
+  // cells_count = 2 ^ log10 ( initial_vertex_count )
+
+  Standard_Real anErrFactorU, anErrFactorV;
+  ComputeErrFactors(myAttribute->GetDefFace (), aSurface, anErrFactorU, anErrFactorV);
+
+  Standard_Integer aCellsCountU, aCellsCountV;
+  if (aType == GeomAbs_Torus)
+  {
+    aCellsCountU = (Standard_Integer)Ceiling(Pow(2, Log10(
+      (myAttribute->GetUMax() - myAttribute->GetUMin()) / myAttribute->GetDeltaX())));
+    aCellsCountV = (Standard_Integer)Ceiling(Pow(2, Log10(
+      (myAttribute->GetVMax() - myAttribute->GetVMin()) / myAttribute->GetDeltaY())));
+  }
+  else if (aType == GeomAbs_Cylinder)
+  {
+    aCellsCountU = (Standard_Integer)Ceiling(Pow(2, Log10(
+      (myAttribute->GetUMax() - myAttribute->GetUMin()) / myAttribute->GetDeltaX() /
+      (myAttribute->GetVMax() - myAttribute->GetVMin()))));
+    aCellsCountV = (Standard_Integer)Ceiling(Pow(2, Log10(
+      (myAttribute->GetVMax() - myAttribute->GetVMin()) / anErrFactorV)));
+  }
+  else
+  {
+    aCellsCountU = (Standard_Integer)Ceiling(Pow(2, Log10(
+      (myAttribute->GetUMax() - myAttribute->GetUMin()) / myAttribute->GetDeltaX() / anErrFactorU)));
+    aCellsCountV = (Standard_Integer)Ceiling(Pow(2, Log10(
+      (myAttribute->GetVMax() - myAttribute->GetVMin()) / myAttribute->GetDeltaY() / anErrFactorV)));
+  }
+
+  AdjustCellsCounts(aSurface, nbVertices, aCellsCountU, aCellsCountV);
+
+  BRepMesh_Delaun trigu(myStructure, tabvert_corr, aCellsCountU, aCellsCountV);
 
   //removed all free edges from triangulation
   const Standard_Integer nbLinks = myStructure->NbLinks();
@@ -376,16 +539,13 @@ void BRepMesh_FastDiscretFace::add(const Handle(BRepMesh_FaceAttribute)& theAttr
     }
   }
 
-  const Handle(BRepAdaptor_HSurface)& gFace = myAttribute->Surface();
-  GeomAbs_SurfaceType thetype = gFace->GetType();
-
   Standard_Boolean rajout = 
-    (thetype == GeomAbs_Sphere || thetype == GeomAbs_Torus);
+    (aType == GeomAbs_Sphere || aType == GeomAbs_Torus);
 
   // Check the necessity to fill the map of parameters
-  const Standard_Boolean useUVParam = (thetype == GeomAbs_Torus         ||
-                                       thetype == GeomAbs_BezierSurface ||
-                                       thetype == GeomAbs_BSplineSurface);
+  const Standard_Boolean useUVParam = (aType == GeomAbs_Torus         ||
+                                       aType == GeomAbs_BezierSurface ||
+                                       aType == GeomAbs_BSplineSurface);
 
   const Standard_Real umax = myAttribute->GetUMax();
   const Standard_Real umin = myAttribute->GetUMin();
@@ -405,12 +565,12 @@ void BRepMesh_FastDiscretFace::add(const Handle(BRepMesh_FaceAttribute)& theAttr
       aDef = control(trigu, Standard_True);
       rajout = (aDef > myAttribute->GetDefFace() || aDef < 0.);
     }
-    if (thetype != GeomAbs_Plane)
+    if (aType != GeomAbs_Plane)
     {
       if (!rajout && useUVParam)
       {
         rajout = (myVParam.Extent() > 2 &&
-          (gFace->IsUClosed() || gFace->IsVClosed()));
+          (aSurface->IsUClosed() || aSurface->IsVClosed()));
       }
 
       if (rajout)

diff --git a/tests/bugs/mesh/bug26889 b/tests/bugs/mesh/bug26889
index fcdea97a1f46f05377bc4cbd390d0633fa5d2ff6..d033335f423c737140c9e2a304556f0527ca749f 100644 (file)
--- a/tests/bugs/mesh/bug26889
+++ b/tests/bugs/mesh/bug26889
@@ -1,5 +1,3 @@
-puts "TODO CR26889 ALL: Error: the mesh takes too long to be built"
-
 puts "========="
 puts "CR26889"
 puts "========="
@@ -19,16 +17,6 @@ dchrono t restart
 incmesh a_1 0.1 1
 dchrono t stop counter incmesh
 
-set info [dlog get]
-dlog reset
-dlog off
-
-regexp {COUNTER incmesh: +([-0-9.+eE]+)} ${info} full Time
-
-if { ${Time} > 30. } {
-  puts "Error: the mesh takes too long to be built"
-}
-
 checktrinfo a_1 -tri -nod
 
 vinit

diff --git a/tests/bugs/mesh/bug29715 b/tests/bugs/mesh/bug29715
new file mode 100644 (file)
index 0000000..f222bba
--- /dev/null
+++ b/tests/bugs/mesh/bug29715
@@ -0,0 +1,15 @@
+puts "======="
+puts "0029715: Mesh - Estimate the grid size of the acceleration structure by the complexity of the face"
+puts "======="
+puts ""
+
+if {[info commands stepread] == ""} {pload XSDRAW}
+
+stepread [locate_data_file bug29715_slow.stp] a *
+renamevar a_1 a
+
+chrono h reset start
+incmesh a 0.0027
+chrono h stop show COUNTER Meshing
+
+tricheck a