0027226: Meshing algorithm fails if edge curves has small defects (bends, loops...)

Take into account tolerance of edge vertices in order to reduce density of discretization points.

src/BRepMesh/BRepMesh_EdgeTessellator.cxx

@@ -30,7 +30,7 @@
 #include <TopTools_ListIteratorOfListOfShape.hxx>
 #include <TopTools_ListOfShape.hxx>
 #include <TColStd_Array1OfReal.hxx>
-
+#include <TopExp.hxx>
 
 IMPLEMENT_STANDARD_RTTIEXT(BRepMesh_EdgeTessellator,BRepMesh_IEdgeTool)
 
@@ -47,6 +47,8 @@ BRepMesh_EdgeTessellator::BRepMesh_EdgeTessellator(
   const Standard_Real                              theMinSize)
   : mySurface(theFaceAttribute->Surface())
 {
+  TopExp::Vertices (theEdge, myFirstVertex, myLastVertex);
+
   Standard_Real aPreciseAngDef = 0.5 * theAngDeflection;
   Standard_Real aPreciseLinDef = 0.5 * theLinDeflection;
   if (theEdge.Orientation() == TopAbs_INTERNAL)
@@ -164,6 +166,20 @@ BRepMesh_EdgeTessellator::BRepMesh_EdgeTessellator(
    myFaceRangeV[1] = mySurface->LastVParameter()  + aDv;
 }
 
+//=======================================================================
+//function : isInToleranceOfVertex
+//purpose  : 
+//=======================================================================
+Standard_Boolean BRepMesh_EdgeTessellator::isInToleranceOfVertex (
+  const gp_Pnt&        thePoint,
+  const TopoDS_Vertex& theVertex)
+{
+  const gp_Pnt        aPoint     = BRep_Tool::Pnt       (theVertex);
+  const Standard_Real aTolerance = BRep_Tool::Tolerance (theVertex);
+
+  return (thePoint.SquareDistance (aPoint) < aTolerance * aTolerance);
+}
+
 //=======================================================================
 //function : Value
 //purpose  : 
@@ -176,6 +192,12 @@ Standard_Boolean BRepMesh_EdgeTessellator::Value(
 {
   myTool->Value(theIndex, mySurface, theParameter, thePoint, theUV);
 
+  if (isInToleranceOfVertex (thePoint, myFirstVertex) ||
+      isInToleranceOfVertex (thePoint, myLastVertex))
+  {
+    return Standard_False;
+  }
+
   // If point coordinates are out of surface range, 
   // it is necessary to re-project point.
   if (mySurface->GetType() != GeomAbs_BSplineSurface &&

src/BRepMesh/BRepMesh_EdgeTessellator.hxx

@@ -24,6 +24,7 @@
 #include <BRepMesh_FaceAttribute.hxx>
 #include <BRepAdaptor_Curve.hxx>
 #include <TopTools_IndexedDataMapOfShapeListOfShape.hxx>
+#include <TopoDS_Vertex.hxx>
 
 class Geom_Surface;
 class Geom2d_Curve;
@@ -76,6 +77,11 @@ private:
                     const Standard_Real         theLast,
                     const Standard_Integer      theNbIter);
 
+  //! Checks whether the given point lies within tolerance of the vertex.
+  Standard_Boolean isInToleranceOfVertex (
+    const gp_Pnt&        thePoint,
+    const TopoDS_Vertex& theVertex);
+
 private:
   NCollection_Handle<BRepMesh_GeomTool> myTool;
   Handle(BRepAdaptor_HSurface)          mySurface;
@@ -85,6 +91,8 @@ private:
   Standard_Real                         myEdgeSqTol;
   Standard_Real                         myFaceRangeU[2];
   Standard_Real                         myFaceRangeV[2];
+  TopoDS_Vertex                         myFirstVertex;
+  TopoDS_Vertex                         myLastVertex;
 };
 
 DEFINE_STANDARD_HANDLE(BRepMesh_EdgeTessellator, BRepMesh_IEdgeTool)

src/GCPnts/FILES

@@ -24,3 +24,7 @@ GCPnts_UniformDeflection.cxx
 GCPnts_UniformDeflection.gxx
 GCPnts_UniformDeflection.hxx
 GCPnts_UniformDeflection.lxx
+GCPnts_DistFunction.hxx
+GCPnts_DistFunction.cxx
+GCPnts_DistFunction2d.hxx
+GCPnts_DistFunction2d.cxx

src/GCPnts/GCPnts_DistFunction.cxx

@@ -0,0 +1,73 @@
+// Copyright (c) 2014-2014 OPEN CASCADE SAS
+//
+// This file is part of Open CASCADE Technology software library.
+//
+// This library is free software; you can redistribute it and/or modify it under
+// the terms of the GNU Lesser General Public License version 2.1 as published
+// by the Free Software Foundation, with special exception defined in the file
+// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
+// distribution for complete text of the license and disclaimer of any warranty.
+//
+// Alternatively, this file may be used under the terms of Open CASCADE
+// commercial license or contractual agreement.
+
+#include <GCPnts_DistFunction.hxx>
+#include <gp_Pnt.hxx>
+
+//=======================================================================
+//function : MaxCurvLinDist
+//purpose  : 
+//=======================================================================
+GCPnts_DistFunction::GCPnts_DistFunction(const Adaptor3d_Curve& theCurve,
+                          const Standard_Real U1, const Standard_Real U2)
+: myCurve(theCurve),
+  myU1(U1), myU2(U2)
+{
+  gp_Pnt P1 = theCurve.Value(U1), P2 = theCurve.Value(U2);
+  myLin = gp_Lin(P1, P2.XYZ() - P1.XYZ());
+}
+//
+//=======================================================================
+//function : Value
+//purpose  : 
+//=======================================================================
+Standard_Boolean GCPnts_DistFunction::Value (const Standard_Real X,
+                                                   Standard_Real& F)
+{
+  if (X < myU1 || X > myU2)
+    return Standard_False;
+  //
+  F = -myLin.SquareDistance(myCurve.Value(X));
+  return Standard_True;
+}
+
+//=======================================================================
+//function : MaxCurvLinDistMV
+//purpose  : 
+//=======================================================================
+
+GCPnts_DistFunctionMV::GCPnts_DistFunctionMV(GCPnts_DistFunction& theCurvLinDist)
+: myMaxCurvLinDist(theCurvLinDist)
+{
+}
+
+//=======================================================================
+//function : Value
+//purpose  : 
+//=======================================================================
+Standard_Boolean GCPnts_DistFunctionMV::Value (const math_Vector& X,
+                                                     Standard_Real& F)
+{
+  Standard_Boolean Ok = myMaxCurvLinDist.Value(X(1), F);
+  return Ok;
+}
+
+//=======================================================================
+//function : NbVariables
+//purpose  : 
+//=======================================================================
+Standard_Integer GCPnts_DistFunctionMV::NbVariables() const
+{
+  return 1;
+}
+

src/GCPnts/GCPnts_DistFunction.hxx

@@ -0,0 +1,67 @@
+// Copyright (c) 2014-2014 OPEN CASCADE SAS
+//
+// This file is part of Open CASCADE Technology software library.
+//
+// This library is free software; you can redistribute it and/or modify it under
+// the terms of the GNU Lesser General Public License version 2.1 as published
+// by the Free Software Foundation, with special exception defined in the file
+// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
+// distribution for complete text of the license and disclaimer of any warranty.
+//
+// Alternatively, this file may be used under the terms of Open CASCADE
+// commercial license or contractual agreement.
+
+#ifndef _GCPnts_DistFunction_HeaderFile
+#define _GCPnts_DistFunction_HeaderFile
+
+#include <gp_Lin.hxx>
+#include <math_Function.hxx>
+#include <math_MultipleVarFunction.hxx>
+#include <Adaptor3d_Curve.hxx>
+
+class gp_Pnt;
+
+//! Class to define function, which calculates square distance between point on curve
+//! C(u), U1 <= u <= U2 and line passing through points C(U1) and C(U2)
+//! This function is used in any minimisation algorithm to define maximal deviation between curve and line,
+//! which required one variable function without derivative (for ex. math_BrentMinimum)
+class GCPnts_DistFunction : public math_Function
+{
+public:
+  Standard_EXPORT GCPnts_DistFunction(const Adaptor3d_Curve& theCurve,
+                                      const Standard_Real U1, const Standard_Real U2);
+  //
+  Standard_EXPORT GCPnts_DistFunction(const GCPnts_DistFunction& theOther);
+
+  Standard_EXPORT virtual Standard_Boolean Value (const Standard_Real X,
+                                                        Standard_Real& F);
+private:
+  GCPnts_DistFunction & operator = (const GCPnts_DistFunction & theOther);
+
+  const Adaptor3d_Curve& myCurve;
+  gp_Lin myLin;
+  Standard_Real myU1;
+  Standard_Real myU2;
+};
+//
+//! The same as class GCPnts_DistFunction, but it can be used in minimization algorithms that
+//! requires multi variable function
+class GCPnts_DistFunctionMV : public math_MultipleVarFunction
+{
+public:
+  Standard_EXPORT GCPnts_DistFunctionMV(GCPnts_DistFunction& theCurvLinDist);
+
+  Standard_EXPORT virtual Standard_Boolean Value (const math_Vector& X,
+                                                        Standard_Real& F);
+
+  Standard_EXPORT virtual Standard_Integer NbVariables() const;
+
+private:
+  GCPnts_DistFunctionMV & operator = (const GCPnts_DistFunctionMV & theOther);
+  GCPnts_DistFunction& myMaxCurvLinDist;
+};
+
+//
+
+
+#endif // _GCPnts_DistFunction_HeaderFile

src/GCPnts/GCPnts_DistFunction2d.cxx

@@ -0,0 +1,76 @@
+// Copyright (c) 2014-2014 OPEN CASCADE SAS
+//
+// This file is part of Open CASCADE Technology software library.
+//
+// This library is free software; you can redistribute it and/or modify it under
+// the terms of the GNU Lesser General Public License version 2.1 as published
+// by the Free Software Foundation, with special exception defined in the file
+// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
+// distribution for complete text of the license and disclaimer of any warranty.
+//
+// Alternatively, this file may be used under the terms of Open CASCADE
+// commercial license or contractual agreement.
+
+#include <GCPnts_DistFunction2d.hxx>
+#include <gp_Pnt2d.hxx>
+
+
+//=======================================================================
+//function : GCPnts_DistFunction2d
+//purpose  : 
+//=======================================================================
+GCPnts_DistFunction2d::GCPnts_DistFunction2d(const Adaptor2d_Curve2d& theCurve,
+                                             const Standard_Real U1, const Standard_Real U2)
+: myCurve(theCurve),
+  myU1(U1), myU2(U2)
+{
+  gp_Pnt2d P2d1 = theCurve.Value(U1), P2d2 = theCurve.Value(U2);
+  myLin = gp_Lin2d(P2d1, P2d2.XY() - P2d1.XY());
+}
+
+//=======================================================================
+//function : Value
+//purpose  : 
+//=======================================================================
+
+Standard_Boolean GCPnts_DistFunction2d::Value (const Standard_Real X,
+                                                     Standard_Real& F)
+{
+  if (X < myU1 || X > myU2)
+    return Standard_False;
+  //
+  gp_Pnt2d aP2d = myCurve.Value(X);
+  F = -myLin.SquareDistance(aP2d);
+  return Standard_True;
+}
+//
+//=======================================================================
+//function : GCPnts_DistFunction2dMV
+//purpose  : 
+//=======================================================================
+GCPnts_DistFunction2dMV::GCPnts_DistFunction2dMV(GCPnts_DistFunction2d& theCurvLinDist)
+: myMaxCurvLinDist(theCurvLinDist)
+{
+}
+
+//=======================================================================
+//function : Value
+//purpose  : 
+//=======================================================================
+Standard_Boolean GCPnts_DistFunction2dMV::Value (const math_Vector& X,
+                                                       Standard_Real& F)
+{
+  Standard_Boolean Ok = myMaxCurvLinDist.Value(X(1), F);
+  return Ok;
+}
+
+
+//=======================================================================
+//function : NbVariables
+//purpose  : 
+//=======================================================================
+Standard_Integer GCPnts_DistFunction2dMV::NbVariables() const
+{
+  return 1;
+}
+

src/GCPnts/GCPnts_DistFunction2d.hxx

@@ -0,0 +1,68 @@
+// Copyright (c) 2014-2014 OPEN CASCADE SAS
+//
+// This file is part of Open CASCADE Technology software library.
+//
+// This library is free software; you can redistribute it and/or modify it under
+// the terms of the GNU Lesser General Public License version 2.1 as published
+// by the Free Software Foundation, with special exception defined in the file
+// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
+// distribution for complete text of the license and disclaimer of any warranty.
+//
+// Alternatively, this file may be used under the terms of Open CASCADE
+// commercial license or contractual agreement.
+
+#ifndef _GCPnts_DistFunction2d_HeaderFile
+#define _GCPnts_DistFunction2d_HeaderFile
+
+#include <gp_Lin2d.hxx>
+#include <math_Function.hxx>
+#include <math_MultipleVarFunction.hxx>
+#include <Adaptor2d_Curve2d.hxx>
+
+class gp_Pnt2d;
+
+//! Class to define function, which calculates square distance between point on curve
+//! C(u), U1 <= u <= U2 and line passing through points C(U1) and C(U2)
+//! This function is used in any minimisation algorithm to define maximal deviation between curve and line,
+//! which required one variable function without derivative (for ex. math_BrentMinimum)
+class GCPnts_DistFunction2d : public math_Function
+{
+public:
+  Standard_EXPORT GCPnts_DistFunction2d(const Adaptor2d_Curve2d& theCurve,
+                                        const Standard_Real U1, const Standard_Real U2);
+  //
+  Standard_EXPORT GCPnts_DistFunction2d(const GCPnts_DistFunction2d& theOther);
+
+  Standard_EXPORT virtual Standard_Boolean Value (const Standard_Real X,
+                                                        Standard_Real& F);
+private:
+  GCPnts_DistFunction2d & operator = (const GCPnts_DistFunction2d & theOther);
+
+  const Adaptor2d_Curve2d& myCurve;
+  gp_Lin2d myLin;
+  Standard_Real myU1;
+  Standard_Real myU2;
+};
+//
+//! The same as class GCPnts_DistFunction2d, 
+//! but it can be used in minimization algorithms that
+//! requires multi variable function
+class GCPnts_DistFunction2dMV : public math_MultipleVarFunction
+{
+public:
+  Standard_EXPORT GCPnts_DistFunction2dMV(GCPnts_DistFunction2d& theCurvLinDist);
+
+  Standard_EXPORT virtual Standard_Boolean Value (const math_Vector& X,
+                                                        Standard_Real& F);
+
+
+  Standard_EXPORT virtual Standard_Integer NbVariables() const;
+
+private:
+  GCPnts_DistFunction2dMV & operator = (const GCPnts_DistFunction2dMV & theOther);
+  GCPnts_DistFunction2d& myMaxCurvLinDist;
+};
+//
+
+
+#endif // _GCPnts_DistFunction2d_HeaderFile

src/GCPnts/GCPnts_TangentialDeflection.cxx

@@ -64,6 +64,22 @@ static void D2 (const Adaptor2d_Curve2d& C, const Standard_Real U,
   VV2.SetCoord (X, Y, 0.0);
 }
 
+static Standard_Real EstimAngl(const gp_Pnt& P1, const gp_Pnt& Pm, const gp_Pnt& P2)
+{
+  gp_Vec V1(P1, Pm), V2(Pm, P2);
+  Standard_Real L = V1.Magnitude() * V2.Magnitude();
+  //
+  if(L > gp::Resolution())
+  {
+    return V1.CrossMagnitude(V2)/L;
+  }
+  else
+  {
+    return 0.;
+  }
+}
+
+
 // Return number of interval of continuity on which theParam is located.
 // Last parameter is used to increase search speed.
 static Standard_Integer getIntervalIdx(const Standard_Real theParam, 
@@ -81,7 +97,7 @@ static Standard_Integer getIntervalIdx(const Standard_Real theParam,
   }
   return anIdx;
 }
-
+//
 //=======================================================================
 //function : CPnts_TangentialDeflection
 //purpose  : 
@@ -161,22 +177,32 @@ Standard_Real GCPnts_TangentialDeflection::ArcAngularStep(
 #include <Geom_BezierCurve.hxx>
 #include <Geom_BSplineCurve.hxx>
 #include <gp_Circ.hxx>
+#include <GCPnts_DistFunction.hxx>
 #define TheCurve Adaptor3d_Curve
 #define Handle_TheBezierCurve   Handle(Geom_BezierCurve)
 #define Handle_TheBSplineCurve  Handle(Geom_BSplineCurve)
+#define TheMaxCurvLinDist GCPnts_DistFunction
+#define TheMaxCurvLinDistMV GCPnts_DistFunctionMV
 #include <GCPnts_TangentialDeflection.gxx>
 #undef Handle_TheBezierCurve
 #undef Handle_TheBSplineCurve
 #undef TheCurve
+#undef TheMaxCurvLinDist
+#undef TheMaxCurvLinDistMV
 
 
 #include <Geom2d_BezierCurve.hxx>
 #include <Geom2d_BSplineCurve.hxx>
 #include <gp_Circ2d.hxx>
+#include <GCPnts_DistFunction2d.hxx>
 #define TheCurve Adaptor2d_Curve2d
 #define Handle_TheBezierCurve   Handle(Geom2d_BezierCurve)
 #define Handle_TheBSplineCurve  Handle(Geom2d_BSplineCurve)
+#define TheMaxCurvLinDist GCPnts_DistFunction2d
+#define TheMaxCurvLinDistMV GCPnts_DistFunction2dMV
 #include <GCPnts_TangentialDeflection.gxx>
 #undef Handle_TheBezierCurve
 #undef Handle_TheBSplineCurve
 #undef TheCurve
+#undef TheMaxCurvLinDist
+#undef TheMaxCurvLinDistMV

src/GCPnts/GCPnts_TangentialDeflection.gxx

@@ -21,6 +21,9 @@
 #include <gp_Pnt.hxx>
 #include <gp_Vec.hxx>
 #include <gp.hxx>
+#include <NCollection_List.hxx>
+#include <math_PSO.hxx>
+#include <math_BrentMinimum.hxx>
 
 #define Us3 0.3333333333333333333333333333
 
@@ -369,6 +372,7 @@ void GCPnts_TangentialDeflection::PerformCurve (const TheCurve& C)
   Standard_Integer aIdx[2] = {Intervs.Lower(), Intervs.Lower()}; // Indexes of intervals of U1 and U2, used to handle non-uniform case.
   Standard_Boolean isNeedToCheck = Standard_False;
   gp_Pnt aPrevPoint = points.Last();
+  Standard_Integer MaxNbCorr = 100;
 
   while (MorePoints) {
     aIdx[0] = getIntervalIdx(U1, Intervs, aIdx[0]);
@@ -383,12 +387,12 @@ void GCPnts_TangentialDeflection::PerformCurve (const TheCurve& C)
     else D0 (C, U2, CurrentPoint);           //Point suivant
 
     Standard_Real Coef = 0.0, ACoef = 0., FCoef = 0.;
-    Standard_Boolean Correction, TooLarge, TooSmall;
+    Standard_Boolean Correction, TooLarge;
     TooLarge   = Standard_False;
     Correction = Standard_True;
-    TooSmall = Standard_False;
+    Standard_Integer nbcorr = 0;
 
-    while (Correction) {                     //Ajustement Du
+    while (Correction && (++nbcorr <= MaxNbCorr)) {                     //Ajustement Du
       if (isNeedToCheck)
       {
         aIdx[1] = getIntervalIdx(U2, Intervs, aIdx[0]);
@@ -447,19 +451,14 @@ void GCPnts_TangentialDeflection::PerformCurve (const TheCurve& C)
           Correction = Standard_False;
         }
         else {
-          if (Coef >= 0.55 || TooLarge) { 
+          if (Coef >= 0.55 || TooLarge || nbcorr == MaxNbCorr) { 
             parameters.Append (U2);
             points    .Append (CurrentPoint);
             aPrevPoint = CurrentPoint;
             Correction = Standard_False;
             isNeedToCheck = Standard_True;
           }
-          else if (TooSmall) {
-            Correction = Standard_False;
-            aPrevPoint = CurrentPoint;
-          }
           else {
-            TooSmall = Standard_True;
             //Standard_Real UUU2 = U2;
             Du += Min((U2-U1)*(1.-Coef), Du*Us3);
 
@@ -478,9 +477,9 @@ void GCPnts_TangentialDeflection::PerformCurve (const TheCurve& C)
             parameters.Append (U1);
             points    .Append (aPrevPoint);
           }
-          U2 = MiddleU;
+          U2 = 0.9*MiddleU + 0.1*U2;
           Du  = U2-U1;
-          CurrentPoint = MiddlePoint;
+          D0 (C, U2, CurrentPoint);
         }
         else {
           Du*=0.9;
@@ -562,4 +561,86 @@ void GCPnts_TangentialDeflection::PerformCurve (const TheCurve& C)
       i++;
     }
   }
+  //Additional check for intervals
+  Standard_Real MinLen2 = myMinLen * myMinLen;
+  Standard_Integer MaxNbp = 10 * Nbp;
+  for(i = 1; i < Nbp; ++i)
+  {
+    U1 = parameters(i);
+    U2 = parameters(i + 1);
+    // Check maximal deflection on interval;
+    Standard_Real dmax = 0.;
+    Standard_Real umax = 0.;
+    Standard_Real amax = 0.;
+    EstimDefl(C, U1, U2, dmax, umax);
+    const gp_Pnt& P1 = points(i);
+    const gp_Pnt& P2 = points(i+1);
+    D0(C, umax, MiddlePoint);
+    amax = EstimAngl(P1, MiddlePoint, P2);
+    if(dmax > curvatureDeflection || amax > AngleMax)
+    {
+      if(umax - U1 > uTol && U2 - umax > uTol)
+      {
+        if (P1.SquareDistance(MiddlePoint) > MinLen2 && P2.SquareDistance(MiddlePoint) > MinLen2)
+        {
+          parameters.InsertAfter(i, umax);
+          points.InsertAfter(i, MiddlePoint);
+          ++Nbp;
+          --i; //To compensate ++i in loop header: i must point to first part of splitted interval
+          if(Nbp > MaxNbp)
+          {
+            break;
+          }
+        }
+      }
+    }
+  }
+  // 
+}
+
+//=======================================================================
+//function : EstimDefl
+//purpose  : Estimation of maximal deflection for interval [U1, U2]
+//           
+//=======================================================================
+void GCPnts_TangentialDeflection::EstimDefl (const TheCurve& C,
+                            const Standard_Real U1, const Standard_Real U2, 
+                            Standard_Real& MaxDefl, Standard_Real& UMax)
+{
+  Standard_Real Du = (lastu - firstu);
+  //
+  TheMaxCurvLinDist aFunc(C, U1, U2);
+ //
+  if(U2 - U1 < 0.01 * Du)
+  {
+    math_BrentMinimum anOptLoc(uTol);
+    anOptLoc.Perform(aFunc, U1, (U1+U2)/2., U2);
+    if(anOptLoc.IsDone())
+    {
+      MaxDefl = Sqrt(-anOptLoc.Minimum());
+      UMax = anOptLoc.Location();
+      return;
+    }
+  }
+  //
+  math_Vector aLowBorder(1,1);
+  math_Vector aUppBorder(1,1);
+  math_Vector aSteps(1,1);
+  //
+  aSteps(1) = Max(0.001 * Du, 10. * uTol);
+  Standard_Integer aNbParticles = Max(8, RealToInt(64 * (U2 - U1) / Du));
+  //
+  aLowBorder(1) = U1;
+  aUppBorder(1) = U2;
+  //
+  //
+  Standard_Real aValue;
+  math_Vector aT(1,1);
+  TheMaxCurvLinDistMV aFuncMV(aFunc);
+
+  math_PSO aFinder(&aFuncMV, aLowBorder, aUppBorder, aSteps, aNbParticles); 
+  aFinder.Perform(aSteps, aValue, aT);
+  MaxDefl = Sqrt(-aValue);
+  UMax = aT(1);
+  //
 }

src/GCPnts/GCPnts_TangentialDeflection.hxx

@@ -27,19 +27,20 @@
 #include <TColStd_SequenceOfReal.hxx>
 #include <Standard_Boolean.hxx>
 #include <gp_Pnt.hxx>
+#include <gp_Lin.hxx>
+#include <math_Function.hxx>
+#include <math_MultipleVarFunction.hxx>
+#include <Adaptor3d_Curve.hxx>
+#include <Adaptor2d_Curve2d.hxx>
 class Standard_ConstructionError;
 class Standard_OutOfRange;
-class Adaptor3d_Curve;
-class Adaptor2d_Curve2d;
-class gp_Pnt;
-
 
 //! Computes a set of  points on a curve from package
 //! Adaptor3d  such  as between  two successive   points
 //! P1(u1)and P2(u2) :
 //!
 //! . ||P1P3^P3P2||/||P1P3||*||P3P2||<AngularDeflection
-//! . ||P1P2^P1P3||/||P1P2||*||P1P3||<CurvatureDeflection
+//! . ||P1P2^P1P3||/||P1P2||<CurvatureDeflection
 //!
 //! where P3 is the point of abscissa ((u1+u2)/2), with
 //! u1 the abscissa of the point P1 and u2 the abscissa
@@ -69,10 +70,12 @@ class gp_Pnt;
 //! U = PointsOnCurve.Parameter (i);
 //! P = PointsOnCurve.Value (i);
 //! }
+
 class GCPnts_TangentialDeflection 
 {
 public:
 
+//
   DEFINE_STANDARD_ALLOC
 
   
@@ -136,6 +139,11 @@ private:
   
   Standard_EXPORT void EvaluateDu (const Adaptor2d_Curve2d& C, const Standard_Real U, gp_Pnt& P, Standard_Real& Du, Standard_Boolean& NotDone) const;
 
+  Standard_EXPORT void EstimDefl (const Adaptor3d_Curve& C, const Standard_Real U1, const Standard_Real U2, 
+                                  Standard_Real& MaxDefl, Standard_Real& UMax);
+  
+  Standard_EXPORT void EstimDefl (const Adaptor2d_Curve2d& C, const Standard_Real U1, const Standard_Real U2, 
+                                  Standard_Real& MaxDefl, Standard_Real& UMax);
 
   Standard_Real angularDeflection;
   Standard_Real curvatureDeflection;

src/GCPnts/GCPnts_UniformDeflection.gxx

@@ -186,7 +186,7 @@ static Standard_Boolean PerformComposite (TColStd_SequenceOfReal& Parameters,
 
     // remove last point to avoid duplication
     Parameters.Remove (Parameters.Length());
-    Points.Remove (Parameters.Length());
+    Points.Remove (Points.Length());
 
     Ua = Ub;
   }

src/GeometryTest/GeometryTest_CurveCommands.cxx

@@ -111,6 +111,7 @@
 #include <GCPnts_QuasiUniformDeflection.hxx>
 #include <GCPnts_UniformDeflection.hxx>
 #include <GCPnts_TangentialDeflection.hxx>
+#include <GCPnts_DistFunction.hxx>
 #include <GeomAPI_ExtremaCurveCurve.hxx>
 #include <gce_MakeLin.hxx>
 #include <TColStd_Array1OfBoolean.hxx>
@@ -831,97 +832,11 @@ static Standard_Integer movelaw (Draw_Interpretor& di, Standard_Integer n, const
 //Static method computing deviation of curve and polyline
 #include <math_PSO.hxx>
 #include <math_PSOParticlesPool.hxx>
-#include <math_MultipleVarFunctionWithHessian.hxx>
-#include <math_NewtonMinimum.hxx>
-
-class aMaxCCDist : public math_MultipleVarFunctionWithHessian
-{
-public:
-  aMaxCCDist(const Handle(Geom_Curve)& theCurve,
-             const Handle(Geom_BSplineCurve)& thePnts)
-: myCurve(theCurve),
-  myPnts(thePnts)
-  {
-  }
-
-  virtual Standard_Boolean Value (const math_Vector& X,
-                                  Standard_Real& F)
-  {
-    if (!CheckInputData(X(1)))
-    {
-      return Standard_False;
-    }
-    F = -myCurve->Value(X(1)).SquareDistance(myPnts->Value(X(1)));
-    return Standard_True;
-  }
-
-  
-  virtual Standard_Boolean Gradient (const math_Vector& X, math_Vector& G)
-  {
-    if (!CheckInputData(X(1)))
-    {
-      return Standard_False;
-    }
-    gp_Pnt aPnt1, aPnt2;
-    gp_Vec aVec1, aVec2;
-    myCurve->D1(X(1), aPnt1, aVec1);
-    myPnts->D1 (X(1), aPnt2, aVec2);
-
-    G(1) = 2 * (aPnt1.X() - aPnt2.X()) * (aVec1.X() - aVec2.X())
-         + 2 * (aPnt1.Y() - aPnt2.Y()) * (aVec1.Y() - aVec2.Y())
-         + 2 * (aPnt1.Z() - aPnt2.Z()) * (aVec1.Z() - aVec2.Z());
-    G(1) *= -1.0; // Maximum search.
-
-    return Standard_True;
-  }
-
-  virtual Standard_Boolean Values (const math_Vector& X, Standard_Real& F, math_Vector& G, math_Matrix& H)
-  {
-    if (Value(X, F) && Gradient(X, G))
-    {
-      gp_Pnt aPnt1, aPnt2;
-      gp_Vec aVec11, aVec12, aVec21, aVec22;
-      myCurve->D2(X(1), aPnt1, aVec11, aVec12);
-      myPnts->D2 (X(1), aPnt2, aVec21, aVec22);
-
-      H(1,1) = 2 * (aVec11.X() - aVec21.X()) * (aVec11.X() - aVec21.X())
-             + 2 * (aVec11.Y() - aVec21.Y()) * (aVec11.Y() - aVec21.Y())
-             + 2 * (aVec11.Z() - aVec21.Z()) * (aVec11.Z() - aVec21.Z())
-             + 2 * (aPnt1.X() - aPnt2.X()) * (aVec12.X() - aVec22.X())
-             + 2 * (aPnt1.Y() - aPnt2.Y()) * (aVec12.Y() - aVec22.Y())
-             + 2 * (aPnt1.Z() - aPnt2.Z()) * (aVec12.Z() - aVec22.Z());
-      H(1,1) *= -1.0; // Maximum search.
-
-      return Standard_True;
-    }
-    return Standard_False;
-  }
-
-  virtual Standard_Boolean Values (const math_Vector& X, Standard_Real& F, math_Vector& G)
-  {
-    return (Value(X, F) && Gradient(X, G));
-  }
-
-  virtual Standard_Integer NbVariables() const
-  {
-    return 1;
-  }
-
-private:
-  aMaxCCDist & operator = (const aMaxCCDist & theOther);
-
-  Standard_Boolean CheckInputData(Standard_Real theParam)
-  {
-    if (theParam < myCurve->FirstParameter() || 
-        theParam > myCurve->LastParameter())
-      return Standard_False;
-    return Standard_True;
-  }
-
-  const Handle(Geom_Curve)& myCurve;
-  const Handle(Geom_BSplineCurve)& myPnts;
-};
+#include <math_MultipleVarFunction.hxx>
+#include <math_BrentMinimum.hxx>
 
+static Standard_Real CompLocalDev(const Handle(Geom_Curve)& theCurve,  
+                                  const Standard_Real u1, const Standard_Real u2);
 
 static void ComputeDeviation(const Handle(Geom_Curve)& theCurve,
                              const Handle(Geom_BSplineCurve)& thePnts,
@@ -939,42 +854,76 @@ static void ComputeDeviation(const Handle(Geom_Curve)& theCurve,
   Standard_Integer nbp = thePnts->NbKnots();
   TColStd_Array1OfReal aKnots(1, nbp);
   thePnts->Knots(aKnots);
-  math_Vector aLowBorder(1,1);
-  math_Vector aUppBorder(1,1);
-  math_Vector aSteps(1,1);
 
   Standard_Integer i;
   for(i = 1; i < nbp; ++i)
   {
-    aLowBorder(1) = aKnots(i);
-    aUppBorder(1) = aKnots(i+1);
-    aSteps(1) =(aUppBorder(1) - aLowBorder(1)) * 0.01; // Run PSO on even distribution with 100 points.
-
-    Standard_Real aValue;
-    math_Vector aT(1,1);
-    aMaxCCDist aFunc(theCurve, thePnts);
-    math_PSO aFinder(&aFunc, aLowBorder, aUppBorder, aSteps); // Choose 32 best points from 100 above.
-    aFinder.Perform(aSteps, aValue, aT);
-    Standard_Real d = 0.;
-
-    math_NewtonMinimum anOptLoc(aFunc);
-    anOptLoc.Perform(aFunc, aT);
-
-    if (anOptLoc.IsDone())
+    Standard_Real u1 = aKnots(i), u2 = aKnots(i+1);
+    Standard_Real d = CompLocalDev(theCurve, u1, u2);
+    if(d > theDmax)
     {
-      d = -anOptLoc.Minimum();
-      if(d > theDmax)
-      {
-        theDmax = d;
-        theUfMax = aLowBorder(1);
-        theUlMax = aUppBorder(1);
-        theImax = i;
-      }
+      theDmax = d;
+      theImax = i;
+      theUfMax = u1;
+      theUlMax = u2;
     }
   }
-  theDmax = Sqrt(theDmax); // Convert to Euclidean distance.
 }
 
+Standard_Real CompLocalDev(const Handle(Geom_Curve)& theCurve,  
+                           const Standard_Real u1, const Standard_Real u2)
+{
+  math_Vector aLowBorder(1,1);
+  math_Vector aUppBorder(1,1);
+  math_Vector aSteps(1,1);
+  GeomAdaptor_Curve TCurve(theCurve);
+  //
+  aLowBorder(1) = u1;
+  aUppBorder(1) = u2;
+  aSteps(1) =(aUppBorder(1) - aLowBorder(1)) * 0.01; // Run PSO on even distribution with 100 points.
+  //
+  GCPnts_DistFunction aFunc1(TCurve,  u1, u2);
+  //
+  Standard_Real aValue;
+  math_Vector aT(1,1);
+  GCPnts_DistFunctionMV aFunc(aFunc1);
+
+  math_PSO aFinder(&aFunc, aLowBorder, aUppBorder, aSteps); // Choose 32 best points from 100 above.
+  aFinder.Perform(aSteps, aValue, aT);
+  Standard_Real d = 0.;
+
+  Standard_Real d1, d2;
+  Standard_Real x1 = Max(u1, aT(1) - aSteps(1));
+  Standard_Boolean Ok = aFunc1.Value(x1, d1);
+  if(!Ok)
+  {
+    return Sqrt(-aValue);
+  }
+  Standard_Real x2 = Min(u2, aT(1) + aSteps(1));
+  Ok = aFunc1.Value(x2, d2);
+  if(!Ok)
+  {
+    return Sqrt(-aValue);
+  }
+  if(!(d1 > aValue && d2 > aValue))
+  {
+    Standard_Real dmin = Min(d1, Min(aValue, d2));
+    return Sqrt(-dmin);
+  }
+
+  math_BrentMinimum anOptLoc(Precision::PConfusion());
+  anOptLoc.Perform(aFunc1, x1, aT(1), x2);
+
+  if (anOptLoc.IsDone())
+  {
+    d = -anOptLoc.Minimum();
+  }
+  else
+  {
+    d = -aValue;
+  }
+  return Sqrt(d);
+}
 
 //=======================================================================
 //function : crvpoints
@@ -1083,6 +1032,7 @@ static Standard_Integer crvtpoints (Draw_Interpretor& di, Standard_Integer n, co
 
   //check deviation
   ComputeDeviation(C,aPnts,dmax,ufmax,ulmax,imax);
+  //
   di << "Max defl: " << dmax << " " << ufmax << " " << ulmax << " " << imax << "\n"; 
 
   return 0;

tests/bugs/iges/buc60820_2

@@ -13,6 +13,6 @@ vdisplay result
 vsetdispmode result 1
 vfit
 
-checktrinfo result -tri 578 -nod 502
+checktrinfo result -tri 618 -nod 533
 
 checkview -display result -2d -path ${imagedir}/${test_image}.png

tests/bugs/iges/buc60823

@@ -14,6 +14,6 @@ vdisplay result
 vsetdispmode result 1
 vfit
 
-checktrinfo result -tri 15571 -nod 9024
+checktrinfo result -tri 18303 -nod 10453
 
 checkview -display result -2d -path ${imagedir}/${test_image}.png

tests/bugs/iges/bug306

@@ -20,7 +20,7 @@ vsetdispmode result 1
 vdisplay result
 vfit
 
-checktrinfo result -tri 19206 -nod 12547
+checktrinfo result -tri 19768 -nod 12899
 
 checkmaxtol result -ref 0.92213088179312575
 checknbshapes result -shell 1

tests/bugs/mesh/bug24127

@@ -14,7 +14,7 @@ incmesh f 1
 
 trinfo f
 
-checktrinfo f -tri 99 -nod 59 -defl 0.59663444648536146 -tol_abs_defl 1.e-3 -tol_rel_defl 0.01
+checktrinfo f -tri 124 -nod 73 -defl 0.39419660561156133 -tol_abs_defl 1.e-3 -tol_rel_defl 0.01
 
 vinit
 vdisplay f

tests/bugs/mesh/bug24938

@@ -1,5 +1,6 @@
-puts "TODO OCC24938 ALL: Error: Number of triangles is equal to 0"
-puts "TODO OCC24938 ALL: Error: Number of nodes is equal to 0"
+#puts "TODO OCC24938 ALL: Error: Number of triangles is equal to 0"
+#puts "TODO OCC24938 ALL: Error: Number of nodes is equal to 0"
+#puts "TODO OCC24938 ALL: Error   : area by triangles differs from the actual area by"
 
 puts "=========="
 puts "OCC24938"

tests/bugs/mesh/bug25519

@@ -15,5 +15,5 @@ fit
 isos a 0
 triangles a
 
-checktrinfo a -tri 2721 -nod 1405 -defl 0.044436924588798624 -tol_rel_defl 0.05 -tol_rel_tri 0.05 -tol_rel_nod 0.05
+checktrinfo a -tri 2792 -nod 1442 -defl 0.041567392896885158 -tol_rel_defl 0.05 -tol_rel_tri 0.05 -tol_rel_nod 0.05
 checkview -screenshot -2d -path ${imagedir}/${test_image}.png

tests/bugs/modalg_2/bug264_7

@@ -15,6 +15,6 @@ vclear
 isos result 0
 triangles result
 
-checktrinfo result -tri 98 -nod 100
+checktrinfo result -tri 106 -nod 108
 checkprops result -s 150.283 
 checkview -display result -3d -path ${imagedir}/${test_image}.png

tests/bugs/modalg_2/bug291

@@ -17,6 +17,6 @@ vfit
 isos result 0
 triangles result
 
-checktrinfo result -tri 1135 -nod 823
+checktrinfo result -tri 1115 -nod 808
 checkprops result -s 376.873 
 checkview -screenshot -3d -path ${imagedir}/${test_image}.png

tests/bugs/modalg_2/bug358

@@ -19,7 +19,7 @@ vdisplay result
 vfit
 vsetdispmode result 1
 
-checktrinfo result -tri 34146 -nod 17507
+checktrinfo result -tri 29010 -nod 14936
 checkprops result -s 24861.2 
 checkshape result
 checkview -display result -2d -path ${imagedir}/${test_image}.png

tests/bugs/modalg_2/bug397

@@ -2,6 +2,9 @@ puts "========================"
 puts " OCC397 "
 puts "========================"
 
+puts "TODO OCC27226 ALL: Colors are not equal in default"
+puts "TODO OCC27226 ALL: Shading is missing in 3D Viewer" 
+
 pload QAcommands
 
 restore [locate_data_file OCC397.brep] a 

tests/bugs/moddata_1/bug15

@@ -14,5 +14,5 @@ vdisplay result
 vsetdispmode result 1
 vfit
 
-checktrinfo result -tri 1009 -nod 593
+checktrinfo result -tri 1043 -nod 612
 checkview -display result -2d -path ${imagedir}/${test_image}.png

tests/bugs/moddata_2/fra62476_2

@@ -13,5 +13,5 @@ tclean result
 incmesh result .1
 triangles result
 
-checktrinfo result -tri 1919 -nod 1008
+checktrinfo result -tri 1807 -nod 952
 checkview -display result -3d -path ${imagedir}/${test_image}.png

tests/bugs/moddata_3/bug25207

@@ -24,17 +24,17 @@ if { ${Nb} != ${expected_Nb} } {
 set tol_abs 1.0e-05
 set tol_rel 0.01
 
-set expected_dmax 0.0013771610718045313
-set expected_ufmax 0.953125
+set expected_dmax 0.00099617706819476581
+set expected_ufmax 0.875
 
 checkreal "dmax" ${dmax} ${expected_dmax} ${tol_abs} ${tol_rel}
 checkreal "ufmax" ${ufmax} ${expected_ufmax} ${tol_abs} ${tol_rel}
 
-set expected_ulmax 0.96875
+set expected_ulmax 0.890625
 if { ${ulmax} != ${expected_ulmax} } {
     puts "Error : bad value of ulmax=${ulmax}"
 }
-set expected_i 73
+set expected_i 68
 if { ${i} != ${expected_i} } {
     puts "Error : bad value of i=${i}"
 }

tests/bugs/moddata_3/bug27108

@@ -0,0 +1,27 @@
+puts "========"
+puts "OCC27108"
+puts "========"
+puts ""
+#################################################################
+# GCPnt_TangentialDeflection does not respect linear deflection
+#################################################################
+
+set bug27108_requested_deflection 0.03
+
+restore [locate_data_file bug27108_Left.brep] a
+
+explode a e
+mkcurve c a_1
+set bug_info [crvtpoints result c $bug27108_requested_deflection 20*pi/180]
+
+smallview
+donly c result
+fit
+
+set bug27108_reached_deflection [lindex $bug_info 6]
+
+if { $bug27108_reached_deflection > $bug27108_requested_deflection} {
+  puts "ERROR: OCC27108 is reproduced. Requested deflection is less than reached."
+}
+
+checkview -screenshot -2d -path ${imagedir}/${test_image}.png

tests/bugs/vis/bug288_5

@@ -13,5 +13,5 @@ isos result 0
 triangles result
 vfit
 
-checktrinfo result -tri 9448 -nod 9080
+checktrinfo result -tri 9202 -nod 8867
 checkview -screenshot -3d -path ${imagedir}/${test_image}.png

tests/bugs/vis/bug364

@@ -18,7 +18,7 @@ vsetdispmode result 1
 vfit
 triangles result
 
-checktrinfo result -tri 92 -nod 92
+checktrinfo result -tri 90 -nod 90
 checkview -screenshot -3d -path ${imagedir}/${test_image}.png
 
 

tests/mesh/data/advanced/A7

@@ -3,5 +3,5 @@ set bug_area "OCC22687"
 set rel_tol 1.3
 if { [string compare $command "shading"] == 0 } {
   set bug_withouttri "OCC22687"
-  set nbwithouttri(ALL) 79
+  set nbwithouttri(ALL) 13
 }

tests/mesh/data/advanced/B6

@@ -3,5 +3,5 @@ if { [string compare $command "incmesh"] == 0 ||
      [string compare $command "mesh"] == 0 ||
      [string compare $command "incmesh_parallel"] == 0 } {
   set bug_area "OCC25519"
-  set rel_tol 1.3485
+  set rel_tol 0.722992148130537 
 }

tests/mesh/data/advanced/B7

@@ -1,5 +1,6 @@
 set TheFileName OCC22302.brep
 set bug_freenodes "OCC22687"
+set bug_area "OCC22687"
 if { [string compare $command "shading"] == 0 } {
   set nbfreenodes(ALL) 4
 } else {

tests/mesh/data/standard/C7

@@ -3,13 +3,16 @@ set bug_freenodes "OCC22687"
 if { [string compare $command "shading"] == 0 } {
    set nbfreenodes(All) 1
 } else {
-   set bug_freelinks "OCC23105"
+   set bug_withouttri "OCC27226"
+###   set bug_freelinks "OCC23105"
 ###   set nbfree(ALL) 4
    if { [string compare $command "mesh"] == 0 } {
 ###     set nbfree(ALL) 8  ### OCC23106
      set nbfree(ALL) 2
+     set nbwithouttri(All) 3
    } else {
      set nbfree(ALL) 2
+     set nbwithouttri(All) 3
    }
    set nbfreenodes(All) 4
 }

tests/mesh/data/standard/O5

@@ -4,6 +4,7 @@ set bug_withouttri "OCC22687"
 set nbwithouttri(ALL) 1
 if { [string compare $command "shading"] == 0 } {
    set rel_tol 0.13
+   set nbwithouttri(ALL) 0
 } else {
    set rel_tol 0.21
 }

tests/mesh/data/standard/U2

@@ -3,12 +3,12 @@ set bug_area "OCC22687"
 set rel_tol 1.9
 set bug_withouttri "OCC22687"
 if { [string compare $command "shading"] == 0 } {
-   puts "TODO OCC23105 ALL: Error: Improvement: The current area difference is"
-   set nbwithouttri(All) 1
+   #puts "TODO OCC23105 ALL: Error: Improvement: The current area difference is"
+   set nbwithouttri(All) 2
    set bug_freenodes "OCC22687"
    set nbfreenodes(All) 38
 } else {
-   set nbwithouttri(All) 1
+   set nbwithouttri(All) 2
    set bug_freenodes "OCC23105"
    set nbfreenodes(ALL) 1
 }

tests/mesh/data/standard/U7

@@ -8,7 +8,7 @@ set bug_freenodes "OCC22687"
 if { [string compare $command "shading"] == 0 } {
 ##set nbt 14
   set nbt 8
-  set nbl 8
+  set nbl 6
   set nbn 83
   set nbwithouttri([checkplatform]) $nbt
   set nbfree([checkplatform]) $nbl

tests/mesh/data/standard/V3

@@ -1 +1,3 @@
 set TheFileName shading_wrongshape_014.brep
+set bug_withouttri "OCC27226"
+set nbwithouttri(ALL) 1

tests/mesh/data/standard/W4

@@ -8,7 +8,7 @@ set bug_freenodes "OCC22687"
 set bug_withouttri "OCC22687"
 if { [string compare $command "shading"] == 0 } {
   set bug_area "OCC22687"
-  set rel_tol 1.3
+  set rel_tol 1.2
   set nbwithouttri(ALL) 6
   set nbfreenodes(ALL) 1
 ##set bug_freelinks "OCC22687"